DATA_SECTION
//Bering Sea Tanner crab model
//********
//********
//to run mcmc 
//scmysr2003bayes -nox -mcmc 1000000 -mcsave 200
// then have to run  scmysr2003bayes -mceval to get output
//whatever is in sd report file will have a distribution and output will go to eval.csv
//for whatever have written to post later in program in the mcmc function part
  number p_const
  int call_no;
  int Nsex;
  int Nmat
  int Nshell
  number spmo // spawning month
  !! spmo=2./12.;
  !!CLASS ofstream post("eval.csv");
  !! p_const = 0.001;
  !! call_no = 0;
  !! Nsex = 2;                                                      // Number of sex
  !! Nmat = 2;                                                      // Number of maturity states
  !! Nshell = 2;                                                    // Number of shell types

  init_int styr                                                      // start year of the model
  init_int endyr                                                    // end year of the model
  !! log_input(styr);
  !! log_input(endyr);

  // Data stuff only from here
  init_int nirec                                                    // number of intial recruitments to estimate 
  init_int nlenm                                                    // number of length bins for males in the model
  init_int nobs_fish                                                // number of years of fishery retained length data
  init_ivector yrs_fish(1,nobs_fish)                                // years when have fishery retained length data
  init_matrix nsamples_fish(1,2,1,nobs_fish)                        // nsamples weight for fish length comps (by year,new/old shell)
  !! log_input(nirec);
  !! log_input(nlenm);
  !! log_input(nobs_fish);
  !! log_input(yrs_fish);
  // Years and sample weights
  // ========================
  init_int nobs_fish_catchf                                          // number of years of directed fishery female and male discard catch data
  init_ivector yrs_fish_catchf(1,nobs_fish_catchf)                   // years which have fishery discard catch data 0 when no fishery

  init_int nobs_fish_discf                                           // number of years of fishery female discard length data
  init_ivector yrs_fish_discf(1,nobs_fish_discf)                     // years which have fishery discard length data
  init_vector nsamples_fish_discf(1,nobs_fish_discf)                 // nsamples weight for fish length comps 
  init_int nobs_fish_discm                                           // number of years of fishery male discard length data
  init_ivector yrs_fish_discm(1,nobs_fish_discm)                     // years which have fishery discard length data
  init_matrix nsamples_fish_discm(1,2,1,nobs_fish_discm)             // nsamples weight for fish length comps (by year and new/old shell)
// snow crab bycatch
  init_int nobs_snowfish_discf                                       // number of years of snow fishery female discard length data
  init_ivector yrs_snowfish_discf(1,nobs_snowfish_discf)             // years which have snow fishery discard length data
  init_vector nsamples_snowfish_discf(1,nobs_snowfish_discf)         // nsamples weight for snow fish length comps 
  init_int nobs_snowfish_discm                                       // number of years of snow fishery male discard length data
  init_ivector yrs_snowfish_discm(1,nobs_snowfish_discm)             // years which have snow fishery discard length data
  init_matrix nsamples_snowfish_discm(1,2,1,nobs_snowfish_discm)     // nsamples weight for snow fish length comps (by year and new/old shell)

  //red king crab bycatch
  init_int nobs_rkfish_discf                                         // number of years of rkc fishery female discard length data
  init_ivector yrs_rkfish_discf(1,nobs_rkfish_discf)                 // years which have rkc fishery discard length data
  init_vector nsamples_rkfish_discf(1,nobs_rkfish_discf)             // nsamples weight for rkc fish length comps 
  init_int nobs_rkfish_discm                                         // number of years of rkc fishery male discard length data
  init_ivector yrs_rkfish_discm(1,nobs_rkfish_discm)                 // years which have rkc fishery discard length data
  init_matrix nsamples_rkfish_discm(1,2,1,nobs_rkfish_discm)         // nsamples weight for rkc fish length comps (by year and new/old shell)
  init_number nobs_discardc                                          // number of years of male directed fishery discard catch data
  init_ivector yrs_discardc(1,nobs_discardc)                         // years which have discard catch data
  //trawl bycatch 
  init_number nobs_trawl_c                                           // number of years of trawl bycatch
  init_ivector yrs_trawl_c(1,nobs_trawl_c)                           // years which have trawl bycatch data
  
  init_number nobs_trawl                                             // number of years of trawl length comps
  init_ivector yrs_trawl(1,nobs_trawl)                               // years which have trawl length data
  init_matrix nsamples_trawl(1,2,1,nobs_trawl)                       // nsamples weight for trawl length comps (by year and ???)
  init_number nobs_srv1                                              // number of years of survey biomass data
  init_ivector yrs_srv1(1,nobs_srv1)                                 // years which have survey biomass estimates
  init_int nobs_srv1_length                                          // number of years of survey male length data
  init_ivector yrs_srv1_length(1,nobs_srv1_length)                   // years which have male length data
  init_4darray nsamples_srv1_length(1,Nmat,1,Nshell,1,Nsex,1,nobs_srv1_length) // number of samples for each length comp by immat,mat,new/old,sex,year

  // for length data
  // first index,1 immat, 2 mature,1 new shell, 2 old shell, then female 1 male 2
  init_5darray obs_p_srv1_lend(1,Nmat,1,Nshell,1,Nsex,1,nobs_srv1_length,1,nlenm)  // immat,mat,new, old survey length data,female,male,year then bin
 LOCAL_CALCS
    log_input(nobs_srv1_length);
    for(int i=1;i<=Nmat;i++)
    {
      for(int j=1;j<=Nshell;j++)
      {
        for(int k=1;k<=Nsex;k++)
        {
          CheckFile<<"#\n"<<"#Mat "<<i<<" Shell "<<j<<" Sex "<<k<<endl;
          for(int l=1;l<=nobs_srv1_length;l++)
          {
            CheckFile<<1973+l<<"\t"<<i<<"\t"<<j<<"\t"<<k<<"\t";
            CheckFile<<obs_p_srv1_lend(i)(j)(k)(l)<<endl;
          }
        }
      }
  }
 END_CALCS
  
  init_3darray obs_p_fish_retd(1,Nshell,1,nobs_fish,1,nlenm)                       // new, old, MALE retained fishery length data
  init_matrix obs_p_fish_discfd(1,nobs_fish_discf,1,nlenm)                         // FEMALE,discard length data
  init_3darray obs_p_fish_discmd(1,Nshell,1,nobs_fish_discm,1,nlenm)               // MALE,discard length data new-old shel
  init_matrix obs_p_snowfish_discf(1,nobs_snowfish_discf,1,nlenm)                         // FEMALE,male snow discard length data
  LOCAL_CALCS
   log_input(nobs_fish);
   log_input(nobs_fish_catchf);
   log_input(nobs_fish_discf);
   log_input(nobs_fish_discm);
   log_input(nobs_snowfish_discf);
   log_input(nobs_snowfish_discm);
   log_input(nobs_rkfish_discf);
   log_input(nobs_rkfish_discm);
   log_input(nobs_trawl_c);
   log_input(nobs_trawl);
   log_input(nobs_srv1);
   log_input(nobs_srv1_length);

   log_input(nsamples_fish);
   log_input(nsamples_fish_discf);
   log_input(nsamples_fish_discm);
   log_input(nsamples_snowfish_discf);
   log_input(nsamples_snowfish_discm);
   log_input(nsamples_rkfish_discm);
   log_input(nsamples_rkfish_discf);
   log_input(nsamples_trawl);
   log_input(yrs_discardc);
   log_input(nobs_srv1_length);
   log_input(yrs_srv1_length);
   log_input(nsamples_srv1_length);
   CheckFile<<" obs survey male mature old shell length "<<endl;
   log_input(obs_p_srv1_lend(2,2,2,nobs_srv1_length));
   CheckFile<<" obs fishery male retained new length "<<endl;
   log_input(obs_p_fish_retd(1,nobs_fish));
   CheckFile<<" obs fishery male retained old length "<<endl;
   log_input(obs_p_fish_retd(2,nobs_fish));
   log_input(obs_p_fish_discfd);
   log_input(obs_p_snowfish_discf);
 END_CALCS
   
  init_3darray obs_p_snowfish_discm(1,2,1,nobs_snowfish_discm,1,nlenm)                         // male new old snow discard length data
  init_matrix obs_p_rkfish_discf(1,nobs_rkfish_discf,1,nlenm)                         // FEMALE red king discard length data
  init_3darray obs_p_rkfish_discm(1,2,1,nobs_rkfish_discm,1,nlenm)                         // male ne oldred king discard length data

  init_3darray obs_p_trawld(1,Nsex,1,nobs_trawl,1,nlenm)                           // male,female trawl discards
  init_vector catch_numbers(1965,endyr)
  init_vector catch_ret(1965,endyr)                                 // retained catch millions of lbs of crab 
  init_matrix catch_odisc(1,2,1,nobs_fish_catchf)                           // estimated directed discard catch numbers female,male
  init_matrix catch_snowodisc(1,2,1,nobs_discardc)                           // estimated snow discard million lbs female,male
  init_matrix catch_rkodisc(1,2,1,nobs_discardc)                           // estimated red king discard millions lbs female,male
  init_vector catch_trawld(1,nobs_trawl_c)                               // trawl bycatch millions lbs sex combined need to apply mort 80%
  init_vector obs_srv1(1,nobs_srv1)                                 // survey numbers (total) in millions of crab
  init_matrix cv_srv1o(1,Nsex,1,nobs_srv1)                          // survey cv
  matrix       cv_srv1(1,2,styr,endyr);
 LOCAL_CALCS
    log_input(nobs_trawl);
    log_input(obs_p_trawld(2));
    log_input(catch_numbers);
    catch_ret       /= 2.2045;                                             //convert to tons
    catch_odisc     /= 2.2045;
    catch_snowodisc /= 2.2045;
    catch_rkodisc   /= 2.2045;
    catch_trawld    /= 2.2045;                                               //convert to 1000s tons
    log_input(catch_ret);
    log_input(catch_odisc);
    log_input(catch_snowodisc);
    log_input(catch_rkodisc);
    log_input(catch_trawld);
  // Survey indices
    obs_srv1 /= 1000.;
    log_input(obs_srv1);
    log_input(cv_srv1o);
 END_CALCS
  
  init_matrix wtf(1,Nmat,1,nlenm)                                   // weight at length juvenile and mature females (from kodiak program)
  init_vector wtm(1,nlenm)                                          // weight at length males (same as used in kodiak)
  
  init_vector maturity_logistic(1,nlenm)                            // logistic maturity probability curve for new shell immature males
  init_matrix maturity_average(1,2,1,nlenm)                         // probability mature for new immature females, avearge proportion mature by length for new males
  init_matrix maturity_old_average(1,2,1,nlenm)                     // average proportion mature by length for old shell females, males
  init_matrix cv_mean_length_obs(1,2,1,2)                           // cv of mean length for female and male min and max age (NOT USED)
  
  init_vector length_bins(1,nlenm)                                  // Midpoints of length bins
  init_vector catch_midpt(styr,endyr)                               // Timing of catches
  init_vector cpue(1969,endyr)                                      // cpue fishery retained males only numbers/potlift
  init_vector catch_ghl(1979,endyr)                                 // historical CPUE data
  init_vector snowcatch(1974,endyr)
  init_vector rkccatch(1969,endyr)
	LOCAL_CALCS
		snowcatch /= 2.2045;
		log_input(snowcatch);
		log_input(length_bins);
		log_input(catch_midpt);
		log_input(cpue);
		log_input(catch_ghl);
		rkccatch /= 2.2045;
		log_input(rkccatch);
		//SJDM trying to extract the relative abundance information that the mdoel is fit to.
		datDump<<"# \n #Relative abundance information"<<endl;
		datDump<<"#Fishery cpue retained males only numbers/potlift \n #Year\t cpue \t cv \t timing"<<endl;
		for(int i=1969;i<=endyr;i++)
		{
			datDump<<i<<"\t"<<cpue(i)<<endl;
		}
		// End of reading normal data file 
		cout << "Completed reading data file" << endl;
		// Open control file....
		ad_comm::change_datafile_name("tc.ctl");
	END_CALCS

  init_int styr_fut                                                 // start year of future projections
  init_int endyr_fut                                                // end year of future projections
  init_number q1                                                    // Q  mult by pop biomass to get survey biomass
  init_vector M_in(1,Nsex)                                          // natural mortality females then males
  init_vector M_matn_in(1,Nsex)                                     // natural mortality mature new shell female/male
  init_vector M_mato_in(1,Nsex)                                     // natural mortality mature old shell female/male
  init_int phase_moltingp                                           // phase to estimate molting prob for mature males
  init_int phase_fishsel                                            // phase to estimate dome shape parameters for fishery selectivities
  init_int growth_switch                                            // switch for which growth function to use
  init_int survsel_som_phase                                            // switch for which survey selectivty to use for 1989 to present - positive estimated negative fixed at somerton and otto
  init_int survsel_phase                                           // switch for fixing all survey sel to somerton and otto - <0 fix, >0 estimate
  int survq_phase
  !! survq_phase = survsel_phase+1;
  init_int phase_fut                                                // phase to do F40% and future projection calculations
  init_vector median_rec(1,2)                                       // median recruitment value to use for last years in model (NOT USED)
  init_int median_rec_yrs                                           // median recruitment fixed for endyr to endyr-median_rec_yrs+1 (NOT USED)
  init_int nsellen                                                  // selectivity is set to the selectivity at nselages-1 after age nselages (NOT USED)
  init_int nsellen_srv1                                             // same as above for survey selectivities
  init_int monot_sel                                                // switch for monotonically increasing selectivities (1 on 0 off) (NOT USED)
  init_int monot_sel_srv1                                           // sames as above for survey (NOT USED)
  init_int phase_logistic_sel                                       // phase to estimate selectivities using logistic function
  init_int phase_selcoffs                                           // phase to estimate smooth selectivities (NOT USED)
  init_vector sel_som(1,5)                                          // parameters for somerton-otto selectivity curve
  
  init_vector wt_like(1,8)                                          // weights for selectivity likelihoods 1 fishery female, 2 survey female, 3 fishery male, 4 survey male
  init_vector like_wght(1,7)                                        // likelihood weights for fishery length data, survey length, age data, catch likelihood, survey biomass likelihood,growth like
  init_number like_wght_mbio                                        // likelihood weight for male biomass fit
  init_number like_wght_rec                                         // ??
  init_number like_wght_recf                                        // ??
  init_number like_wght_sexr                                        // ??
  init_number like_wght_sel50                                       // ??
  init_number like_wght_fph1                                        // ??
  init_number like_wght_fph2                                        // ??
  init_number like_wght_fdev                                        // ??
    //weights for constraint on monotonicity 1 fishery female, 2 fishery male, 3 survey female, 4 survey male

  init_number m_disc                                                // fraction of pot discards that die (.5)
  init_number m_trawl                                               // fraction of trawl discards that die(.8)
  
  init_number linff_obs                                             // Growth pars #1 (NOT USED
  init_number sd_linff                                              // Growth pars #1
  init_number linfm_obs                                             // Growth pars #1
  init_number sd_linfm                                              // Growth pars #1
  init_number growthkf_obs                                          // Growth pars #1
  init_number sd_growthkf                                           // Growth pars #1
  init_number growthkm_obs                                          // Growth pars #1
  init_number sd_growthkm                                           // Growth pars #1
  
  init_number af1_obs                                                // A-female
  init_number af2_obs                                                // A-female
  init_number sd_af                                                 // SD (A)
  init_number am1_obs                                                // A-male
  init_number am2_obs                                                // A-male
  init_number sd_am                                                 // SD (A)
  init_number bf1_obs                                                // B-female
  init_number bf2_obs                                                // B-female
  init_number sd_bf                                                 // SD (B) 
  init_number bm1_obs                                                // B-male
  init_number bm2_obs                                                // B-male
  init_number sd_bm                                                 // SD (B)
  
  init_number var_rec_obs                                           // (NOT USED)
  init_number sd_var_rec                                            // (NOT USED)
  init_number var_last_obs                                          // (NOT USED)
  init_number sd_var_last                                           // (NOT USED)
  init_number mate_ratio                                            // mating ratio (USED)
  init_number fraction_new_error                                    // accounts for shell error
  init_number maturity_switch                                       // Set > 0 for logistic maturity instead of fractions by year (males)
  init_int nages                                                    // number of ages to track for mature old shell 
  init_number wght_total_catch                                      // weight for total catch biomass (WHY here!)
  init_number wght_female_potcatch                                  // weight for female pot bycatch
  init_number cpue_cv                                               // cv for fit to fishery pot cpue
  init_number  wt_lmlike
  init_int     q_prior_switch
  init_int     killem_switch
  init_int     n_qs                                                  // number of catchability periods
  init_ivector yrs_q(1,n_qs)                                     // years in which catchability changed
  init_number  end_of_ctl;
  !!cout<<"end of data reading"<<endl;
  int styr_rec;                                                     // rest are working variables 

 LOCAL_CALCS
  log_input(styr_fut);
  log_input(endyr_fut);
  log_input(q1);
  log_input(M_in(1,Nsex));
  log_input(M_matn_in(1,Nsex));
  log_input(M_mato_in(1,Nsex));
  log_input(phase_moltingp);
  log_input(phase_fishsel);
  log_input(growth_switch );
  log_input(survsel_som_phase  );
  log_input(survsel_phase);
  log_input(phase_fut);
  log_input(median_rec(1,2));
  log_input(median_rec_yrs);
  log_input(nsellen);
  log_input(nsellen_srv1);
  log_input(monot_sel);
  log_input(monot_sel_srv1);
  log_input(phase_logistic_sel);
  log_input(phase_selcoffs);
  log_input(sel_som(1,5));
  log_input(wt_like(1,8));
  log_input(like_wght(1,7));
  log_input(like_wght_mbio);
  log_input(like_wght_rec);
  log_input(like_wght_recf);
  log_input(like_wght_sexr);
  log_input(like_wght_sel50);
  log_input(like_wght_fph1);
  log_input(like_wght_fph2);
  log_input(like_wght_fdev);
  log_input(cpue_cv);
  log_input(wt_lmlike);
  log_input(end_of_ctl);
  CheckFile <<" nages "<<nages<<endl;
  // Compute CV's to use (with likelihoods and annual varying CVs) REMOVED VARIANCE INFLATION FACTORS...) Jim
  for(int i=1;i<=nobs_srv1;i++)
  {
    cv_srv1(1,yrs_srv1(i)) = cv_srv1o(1,i);
    cv_srv1(2,yrs_srv1(i)) = cv_srv1o(2,i);
  }

   styr_rec = styr-nirec;                                           // year to start estimating recruits to get initial age comp
   if(nsellen>nlenm) nsellen=nlenm;                                 // make sure nselages not greater than nages
   if(nsellen_srv1>nlenm) nsellen_srv1=nlenm;                       // same as above for survey
   obs_srv1=obs_srv1*1000000;                                       // survey numbers read in are millions of crab
   wtf=wtf*0.001;                                                   // change weights to tons
   wtm=wtm*0.001;                                                    // change weights to tons
   catch_odisc = m_disc*catch_odisc;                                //apply discard mortality to catches
   catch_snowodisc = m_disc*catch_snowodisc;
   catch_rkodisc = m_disc*catch_rkodisc;
   catch_trawld = m_trawl*catch_trawld;
 END_CALCS

  vector tmps(1,nlenm)                                              // Temporrary variable
  vector sumsrv(1,nobs_srv1_length)                                 // Total survey numbers
  
  vector catch_tot(styr,endyr);                                     // Total catches
  matrix catch_disc(1,2,styr,endyr)                                 // Discard catch
  3darray obs_p_fish_ret(1,2,1,nobs_fish,1,nlenm)                   // length-frequency of retained catch?
  3darray obs_p_fish_tot(1,2,1,nobs_fish,1,nlenm)                   // length-frequency of total catch?
  3darray obs_p_fish_discm(1,2,1,nobs_fish_discm,1,nlenm)           // males discards
  matrix obs_p_fish_discf(1,nobs_fish_discf,1,nlenm)                // female discards
  3darray obs_p_trawl(1,2,1,nobs_trawl,1,nlenm)                     // trawl discards
  5darray obs_p_srv1_len(1,2,1,2,1,2,1,nobs_srv1_length,1,nlenm)    // Length-frequency (survey)
  5darray obs_p_srv1_len1(1,2,1,2,1,2,1,nobs_srv1_length,1,nlenm) 
  matrix obs_srv1_spbiom(1,2,styr,endyr)                            // Survey biomass (by sex)
  3darray obs_srv1_spnum(1,2,1,2,styr,endyr)                        // Survey number (by sex)
  3darray obs_p_snow(1,2,1,nobs_snowfish_discf,1,nlenm)
  3darray obs_p_rk(1,2,1,nobs_rkfish_discf,1,nlenm)
  vector obs_lmales(1,nobs_srv1_length)                             // Large males in survey
  vector obs_lmales_bio(1,nobs_srv1_length)                         // Male biomass
  
  3darray obs_srv1_num(1,2,styr,endyr,1,nlenm)                      // Survey numbers
  vector obs_srv1t(styr,endyr)                                      // ???
  matrix obs_srv1_bioms(1,2,styr,endyr)                             // Survey biomass by sex
  vector obs_srv1_biom(styr,endyr)                                  // Total survey biomass

  vector avgpf(1,nlenm)
  vector avgpm(1,nlenm)
  vector avgp(1,nlenm)
  vector obs_catcht_biom(styr,endyr)                                // Total catch  
  vector obs_catchdm_biom(styr,endyr)                               // Male discards
  vector obs_catchdf_biom(styr,endyr)                               // Female discards
  vector obs_catchtot_biom(styr,endyr)                              // Total catch (discard and retained_
  number avgwt2
  number avgwtall
  vector avgwt(styr,endyr)

// =======================================================================

INITIALIZATION_SECTION
  mean_log_rec1 11.4
  af1 0.5656
  bf1 0.913266
  am1 0.437941
  bm1 0.9487
  log_avg_fmort -0.7
//  log_avg_fmortdf -1.0
  log_avg_fmortt  -4.0
  log_avg_fmortd_snow -3.0
  log_avg_sel50_mn  4.87  //this is 130.3 mm
  fmort_dev 0.00001
  matestf -1.0
  matestm -1.0
//  fish_disc_slope_tf 0.05
//  fish_disc_sel50_tf 85.0
//  fish_disc_slope_tm 0.07
//  fish_disc_sel50_tm 65.0
  
//  moltp_af 2.0
//  moltp_bf 150.
//  moltp_am 0.02
//  moltp_bm 300.
//  moltp_ammat 0.05
//  moltp_bmmat 105.0
    
//  af 15.75
//  bf 1.01
//  am2 15.75
//  bm2 1.07
  
//  srv1_slope 0.07
//  srv1_sel50 60.0
//  srv1_sel95 100
//  srv1_sel50 60 
//  srv2_sel95 100
//  srv2_sel50  60
//  srv3_sel95 100
//  srv3_sel50  60
//  fish_fit_sel50_mn 95.1
//  log_sel50_dev_mn  0.0000
  cpueq 0.001

// =======================================================================

PARAMETER_SECTION
  init_bounded_number af1(0.4,0.77,-8)                       // Female growth-increment
  init_bounded_number bf1(0.6,1.2,-8)                       // Female growth-increment
  init_bounded_number am1(0.3,0.6,-8)                       // Male growth-increment
  init_bounded_number bm1(0.7,1.2,-8)                       // Male growth-increment
  init_bounded_vector growth_beta(1,Nsex,0.75000,0.75001,-2)        // Growth beta
  init_bounded_vector matestf(1,16,-15.0,0.0,5)
  init_bounded_vector matestm(1,nlenm,-15.0,0.0,5)
  init_bounded_number Mmult_imat(0.2,2.0,7)                         // natural mortality females and males
  init_bounded_number Mmultm(0.1,1.9,7)                              // natural mortality mature new and old shell male
  init_bounded_number Mmultf(0.1,1.9,7)                             // natural mortality mature new and old shell female
  init_bounded_vector mat_big(1,2,0.1,10.0,8)                  //mult. on 1983 M for mature males and females                     
  init_bounded_number alpha1_rec(11.49,11.51,-8)               // Parameters related to fraction recruiting
  init_bounded_number beta_rec(3.99,4.01,-8)                   // Parameters related to fraction recruiting
  
  init_bounded_number moltp_af(0.04,3.0,-6)                    // paramters for logistic function molting
  init_bounded_number moltp_bf(130.,300.,-6)                   // female
  init_bounded_number moltp_am(0.04,3.0,-5)                    // paramters for logistic function molting
  init_bounded_number moltp_bm(130.0,300.0,-5)                 // immature males
  init_bounded_number moltp_ammat(.0025,3.0,phase_moltingp)    // logistic molting prob for mature males
  init_bounded_number moltp_bmmat(1,120,phase_moltingp)        // logistic molting prob for mature males

  init_number mean_log_rec1(1)                                 // Mean recruitment
  init_number mean_log_rec1_early(1)                                 // Mean recruitment
  vector mean_log_rec(1,2)  
  init_bounded_dev_vector rec_devf(1974,endyr-1,-15,15,1)      // Deviations about mean recruitment
  init_bounded_dev_vector rec_devm(1974,endyr-1,0,0,-3)
  init_bounded_dev_vector rec_devf_early(styr,1973,-15,15,1)      // Deviations about mean recruitment
  //directed fishery closed 10 years                                 // fishing mortality (directed fishery)
  init_number log_avg_fmort(1)
  init_bounded_vector fmort_dev(1966,endyr-12,-15,15,2)
  init_number log_avg_fmortt(1)                                   // fishing mortality (trawl)
  init_bounded_vector fmortt_dev(1973,endyr,-15,15,3)
  init_number log_avg_fmortd_snow(3)                               // fishing mortality snow crab fishery discards
  init_bounded_vector fmortd_snow_dev(1992,endyr,-15,15,4)
  init_bounded_number log_avg_fmortd_rk(-5.25,-5.25,-4)                               // fishing mortality red king crab fishery discards
  init_bounded_vector fmortd_rk_dev(1,nobs_discardc-1,-15,15,-5)
  // Selectivity pattern for males (directed fishery)
  // Set -phase so not estimated if using @3 selectivity periods
  init_bounded_number fish_slope_mn(0.1,0.4,-phase_logistic_sel)      
  init_bounded_number log_avg_sel50_mn(4,5.0,-phase_logistic_sel)
  init_bounded_vector log_sel50_dev_mn(1,nobs_fish,-5,5,-phase_logistic_sel)
  vector fish_sel50_mn(styr,endyr)
  
// Retention function
// init_bounded_number fish_fit_slope_mn(.250,1.001,phase_logistic_sel)
// init_bounded_number fish_fit_sel50_mn(85.0,160.0,phase_logistic_sel)
// 1981 - 1992
  init_bounded_number fish_fit_slope_mn1(.250,1.001,phase_logistic_sel)
  init_bounded_number fish_fit_sel50_mn1(85.0,160.0,phase_logistic_sel)
// 2005-endyr  
  init_bounded_number fish_fit_slope_mn2(.250,2.001,phase_logistic_sel)
  init_bounded_number fish_fit_sel50_mn2(85.0,160.0,phase_logistic_sel)

  // Directed fishery selectivity pattern for period-1: 1993-1996
  init_bounded_number fish_slope_1(0.05,0.75,phase_logistic_sel)      
  init_bounded_number fish_sel50_1(50,170,-phase_logistic_sel)
  
  // Directed fishery selectivity pattern for period-2: 2005-p
//  init_bounded_number fish_slope_2(0.1,0.5,phase_logistic_sel)      
//  init_bounded_number fish_sel50_2(50,170,phase_logistic_sel)
  
  // Directed fishery selectivity pattern for period-3: 2005-P
//  init_bounded_number fish_slope_3(0.01,0.5,-phase_logistic_sel)      
//  init_bounded_number fish_sel50_3(50,170,-phase_logistic_sel)

    // Directed fishery selectivity pattern changing by year for period-3: 2005-P
   init_bounded_number fish_slope_yr_3(0.1,0.4,phase_logistic_sel)      
  init_bounded_number log_avg_sel50_3(4,5.0,phase_logistic_sel)
  init_bounded_dev_vector log_sel50_dev_3(1,11,-.5,.5,phase_logistic_sel)
  matrix sel_byyr_3(2005,endyr,1,nlenm)
  
  // for a dome-shaped selex pattern
  init_bounded_number fish_slope_mn2(.01,2.0,phase_fishsel)
  init_bounded_number fish_sel50_mn2(100.0,160.0,phase_fishsel)

  // Female discards
  init_bounded_number fish_disc_slope_f(.1,.4,phase_logistic_sel)
  init_bounded_number fish_disc_sel50_f(80.0,150.0,phase_logistic_sel)

   // snow fishery female discards
//  init_bounded_number snowfish_disc_slope_f(.05,.5,phase_logistic_sel+1)         //add 1 to phase
//  init_bounded_number snowfish_disc_sel50_f(50.0,100.0,phase_logistic_sel+1)

   // snow fishery female discards for period-1: 1989-1996
  init_bounded_number snowfish_disc_slope_f_1(.05,.5,phase_logistic_sel+1)         //add 1 to phase
  init_bounded_number snowfish_disc_sel50_f_1(50.0,150.0,phase_logistic_sel+1)

   // snow fishery female discards for period-2: 1997-2004
  init_bounded_number snowfish_disc_slope_f_2(.05,.5,phase_logistic_sel+1)         //add 1 to phase
  init_bounded_number snowfish_disc_sel50_f_2(50.0,120.0,phase_logistic_sel+1)

   // snow fishery female discards for period-3: 2005-P
  init_bounded_number snowfish_disc_slope_f_3(.05,.5,phase_logistic_sel+1)         //add 1 to phase
  init_bounded_number snowfish_disc_sel50_f_3(50.0,120.0,phase_logistic_sel+1)
  
  // snow fishery male discards
//  init_bounded_number snowfish_disc_slope_m(.1,.5,phase_logistic_sel+1)          //1 to phase
//  init_bounded_number snowfish_disc_sel50_m(60.0,150.0,phase_logistic_sel+1)
//  init_bounded_number snowfish_disc_slope_m2(.1,.5,phase_logistic_sel+1)
//  init_bounded_number snowfish_disc_sel50_m2(40.0,200.0,phase_logistic_sel+1)
  
// snow fishery male discards for period-1: 1989-1996
  init_bounded_number snowfish_disc_slope_m_1(.1,.5,phase_logistic_sel+1)          //1 to phase
  init_bounded_number snowfish_disc_sel50_m_1(60.0,150.0,phase_logistic_sel+1)
  init_bounded_number snowfish_disc_slope_m2_1(.1,.5,phase_logistic_sel+1)
  init_bounded_number snowfish_disc_sel50_m2_1(40.0,200.0,phase_logistic_sel+1)
  
// snow fishery male discards for period-2: 1997-2004
  init_bounded_number snowfish_disc_slope_m_2(.1,.5,phase_logistic_sel+1)          //1 to phase
  init_bounded_number snowfish_disc_sel50_m_2(60.0,150.0,phase_logistic_sel+1)
  init_bounded_number snowfish_disc_slope_m2_2(.1,.5,phase_logistic_sel+1)
  init_bounded_number snowfish_disc_sel50_m2_2(40.0,200.0,phase_logistic_sel+1)
  
// snow fishery male discards for period-3: 2005-P
  init_bounded_number snowfish_disc_slope_m_3(.1,.5,phase_logistic_sel+1)          //1 to phase
  init_bounded_number snowfish_disc_sel50_m_3(60.0,150.0,phase_logistic_sel+1)
  init_bounded_number snowfish_disc_slope_m2_3(.1,.5,phase_logistic_sel+1)
  init_bounded_number snowfish_disc_sel50_m2_3(40.0,200.0,phase_logistic_sel+1)
  
    // red king fishery female discards
//  init_bounded_number rkfish_disc_slope_f(.05,.5,-phase_logistic_sel)   //add 2 to phase
//  init_bounded_number rkfish_disc_sel50_f(75.0,115.0,-phase_logistic_sel)   //add 2 to phase
  init_bounded_number rkfish_disc_slope_f1(.05,.5,phase_logistic_sel)   //add 2 to phase
  init_bounded_number rkfish_disc_sel50_f1(50.0,150.0,phase_logistic_sel)   //add 2 to phase
  init_bounded_number rkfish_disc_slope_f2(.05,.5,phase_logistic_sel)   //add 2 to phase
  init_bounded_number rkfish_disc_sel50_f2(50.0,150.0,phase_logistic_sel)   //add 2 to phase
  init_bounded_number rkfish_disc_slope_f3(.05,.5,phase_logistic_sel)   //add 2 to phase
  init_bounded_number rkfish_disc_sel50_f3(50.0,170.0,phase_logistic_sel)   //add 2 to phase
  
    // red king fishery male discards
//  init_bounded_number rkfish_disc_slope_m(.05,.30,-phase_logistic_sel)          //add 2 to phase
//  init_bounded_number rkfish_disc_sel50_m(95.0,125.0,-phase_logistic_sel)
  init_bounded_number rkfish_disc_slope_m1(.01,.50,phase_logistic_sel)          //add 2 to phase
  init_bounded_number rkfish_disc_sel50_m1(95.0,150.0,phase_logistic_sel)
  init_bounded_number rkfish_disc_slope_m2(.01,.50,phase_logistic_sel)          //add 2 to phase
  init_bounded_number rkfish_disc_sel50_m2(95.0,150.0,phase_logistic_sel)
  init_bounded_number rkfish_disc_slope_m3(.01,.50,phase_logistic_sel)          //add 2 to phase
  init_bounded_number rkfish_disc_sel50_m3(95.0,150.0,phase_logistic_sel)

  // Trawl fishery selectivity female, 1973-1987
  init_bounded_number fish_disc_slope_tf1(0.01,0.5,phase_logistic_sel)
  init_bounded_number fish_disc_sel50_tf1(40.0,125.01,phase_logistic_sel)
  // Trawl fishery selectivity female, 1988-1996
  init_bounded_number fish_disc_slope_tf2(0.005,0.5,phase_logistic_sel)
  init_bounded_number fish_disc_sel50_tf2(40.0,250.01,phase_logistic_sel) 
  // Trawl fishery selectivity female, 1997-P
  init_bounded_number fish_disc_slope_tf3(0.01,0.5,phase_logistic_sel)
  init_bounded_number fish_disc_sel50_tf3(40.0,150.01,phase_logistic_sel)
  // Trawl fishery selectivity male, 1973-1987
  init_bounded_number fish_disc_slope_tm1(0.01,0.5,phase_logistic_sel)
  init_bounded_number fish_disc_sel50_tm1(40.0,120.01,phase_logistic_sel)
  // Trawl fishery selectivity male, 1988-1996
  init_bounded_number fish_disc_slope_tm2(0.01,0.5,phase_logistic_sel)
  init_bounded_number fish_disc_sel50_tm2(40.0,120.01,phase_logistic_sel)
  // Trawl fishery selectivity male, 1997-P
  init_bounded_number fish_disc_slope_tm3(0.01,0.5,phase_logistic_sel)
  init_bounded_number fish_disc_sel50_tm3(40.0,120.01,phase_logistic_sel)

  init_bounded_number srv_sel95_m(70.0,200.,survsel_phase)
  init_bounded_number srv_sel50_m(0.0,90.,survsel_phase)
  init_bounded_number srv_sel95_f(70.0,200.,survsel_phase)
  init_bounded_number srv_sel50_f(0.0,90.,survsel_phase)
  init_vector         srv_sel95_m_dev(1,n_qs-1,survsel_phase)
  init_vector         srv_sel50_m_dev(1,n_qs-1,survsel_phase)
  init_vector         srv_sel95_f_dev(1,n_qs-1,survsel_phase)
  init_vector         srv_sel50_f_dev(1,n_qs-1,survsel_phase)
    
  init_bounded_number srv_q_f_par(0.0,1.2,survsel_phase)
  init_bounded_number srv_q_m_par(0.0,1.2,survsel_phase)
  init_vector         srv_q_f_dev(1,n_qs-1,survq_phase)
  init_vector         srv_q_m_dev(1,n_qs-1,survq_phase)

  // pot fishery cpue q
  init_bounded_number cpueq(0.00055,0.00055,-5)
  
  // proportion new shell in recruits
  init_bounded_number proprecn(1.0,1.0,-2)
  ////end of estimated parameters///////////////

  3darray sel(1,2,styr,endyr,1,nlenm)                               // Directed selectivity (males)
  vector self_use(1,nlenm)           // Temporary vector for @3 selectivity periods
  vector sel_1(1,nlenm)            // Directed fishery selectivity (males): period-1
  vector sel_2(1,nlenm)            // Directed fishery selectivity (males): period-2
  vector sel_3(1,nlenm)            // Directed fishery selectivity (males): period-3
  vector sel_ret1(1,nlenm)                 // Retention function (fraction crabxlength retained)
  vector sel_discf(1,nlenm)                                         // Directed selectivity (females)
  3darray sel_fit(1,2,styr,endyr,1,nlenm)                           // Net selectivity (gear + retained)
  3darray sel_ret(1,2,styr,endyr,1,nlenm)                           // Retension curve
//  matrix sel_trawl(1,2,1,nlenm)                                   // Trawl selectivity
  3darray sel_trawl(1,3,1,2,1,nlenm)        // Trawl selectivity (time period x sex x length bin)
  matrix sel_trawl_use(1,2,1,nlenm)
//  matrix sel_disc_snow(1,2,1,nlenm)
  3darray sel_disc_snow(1,3,1,2,1,nlenm)       // 3D array to accomated 3 selectivity periods
  matrix sel_disc_snow_use(1,2,1,nlenm)

//    matrix sel_disc_rkc(1,2,1,nlenm)
  3darray sel_disc_rkc(1,3,1,2,1,nlenm)              // 3D array to accomated 3 selectivity periods
  matrix sel_disc_rkc_use(1,2,1,nlenm)
  
  matrix sel_srv_f(styr,endyr,1,nlenm)                            // Survey selectivity 1
  matrix sel_srv_m(styr,endyr,1,nlenm)                            // Survey selectivity 1
  vector srv_q_f(styr,endyr)                                      // Survey q 
  vector srv_q_m(styr,endyr)                                      // Survey q 


  vector popn(styr,endyr)                                           // Total population numbers (output)
  number popn_snowm                                           // Total population numbers (output)
  number popn_snowf                                           // Total population numbers (output)
  number popn_rkm                                           // Total population numbers (output)
  number popn_rkf                                           // Total population numbers (output)

  vector M(1,2)
  vector M_matn(1,2)
  vector M_mato(1,2)
  
    vector pred_bio(styr,endyr)                                       // Predicted biomass (determines depletion)
  
  matrix totn_srv1(1,2,styr,endyr)                                  // total survey abundance
  vector fspbio(styr,endyr)                                         // Predicted female survey biomass
  vector mspbio(styr,endyr)                                         // Predicted male survey biomass
  vector legal_males(styr,endyr)                                    // Legal number of males
  vector legal_males_bio(styr,endyr)                                // Legal biomass (output)
  vector legal_srv_males(styr,endyr)                                // Survey-selected males (output)
  vector legal_srv_males_n(styr,endyr)                              // Survey-selected males (output)
  vector legal_srv_males_o(styr,endyr)                              // Survey-selected males (output)
  vector legal_srv_males_bio(styr,endyr)                            // Survey-selected males (output)
  matrix biom_tmp(1,2,styr,endyr);                                  // Predicted survey indices
  vector fspbio_srv1(styr,endyr)                                    // Female survey biomass
  vector mspbio_srv1(styr,endyr)                                    // Male survey biomass
  matrix fspbio_srv1_num(1,2,styr,endyr)                            // Survey biomass (females) - output
  matrix mspbio_srv1_num(1,2,styr,endyr)                            // Survey biomass (males) - output
  matrix pred_srv1_bioms(1,2,styr,endyr)                            // Survey biomass (an output)
  3darray pred_srv1(1,2,styr,endyr,1,nlenm)                         // Predicted survey - output
  4darray pred_p_srv1_len_new(1,2,1,2,styr,endyr,1,nlenm)           // Predicted new shell length-frequency
  4darray pred_p_srv1_len_old(1,2,1,2,styr,endyr,1,nlenm)           // Predicted old shell length-frequency
  
  3darray pred_p_fish(1,2,styr,endyr,1,nlenm)                       // Predicted proportion (total catch)
  3darray pred_p_fish_fit(1,2,styr,endyr,1,nlenm)                   // Predicted retained catch proportions
  matrix  pred_p_fish_discf(styr,endyr,1,nlenm)                     // Predicted female discard proprtions
  3darray pred_p_trawl(1,2,styr,endyr,1,nlenm)                      // Predicted tarwl proportions
  3darray pred_p_snow(1,2,styr,endyr,1,nlenm)                      // Predicted snow crab fishery  proportions
  3darray pred_p_rk(1,2,styr,endyr,1,nlenm)                      // Predicted red king crab proportions
  
  vector catch_trawl(styr,endyr)
  matrix catch_lmale_new(styr,endyr,1,nlenm)                        // Predicted catch (males, total, new shell)
  matrix catch_lmale_old(styr,endyr,1,nlenm)                        // Predicted catch (males, total, old shell)
  matrix catch_male_ret_new(styr,endyr,1,nlenm)                     // Predicted catch (males, retained, new shell)
  matrix catch_male_ret_old(styr,endyr,1,nlenm)                     // Predicted catch (males, retained, old shell)  
  vector pred_catch(styr,endyr)                                     // Total catch (males, directed)
  vector pred_catch_ret(styr,endyr)                                 // Retained catch (males, directed)
  matrix pred_catch_disc(1,2,styr,endyr)                            // Discard catch (in mass by sex)
  vector pred_catch_trawl(styr,endyr)                               // Trawl catch (in mass and by sex)
  vector pred_catch_snowd(styr,endyr)
  vector pred_catch_female_snowd(styr,endyr)
  vector pred_catch_female_d(styr,endyr)
  vector pred_catch_rkd(styr,endyr)
  vector pred_catch_female_rkd(styr,endyr)

  matrix catch_female_d_new(styr,endyr,1,nlenm)
  matrix catch_female_d_old(styr,endyr,1,nlenm)
  matrix catch_female_d(styr,endyr,1,nlenm)

  matrix catch_male_snowd(styr,endyr,1,nlenm)
  matrix catch_female_snowd(styr,endyr,1,nlenm)
  matrix catch_male_snowd_new(styr,endyr,1,nlenm)
  matrix catch_male_snowd_old(styr,endyr,1,nlenm)
  matrix catch_male_rkd_new(styr,endyr,1,nlenm)
  matrix catch_male_rkd_old(styr,endyr,1,nlenm)
  matrix catch_male_rkd(styr,endyr,1,nlenm)
  matrix catch_female_rkd(styr,endyr,1,nlenm)
  matrix catch_female_snowd_new(styr,endyr,1,nlenm)
  matrix catch_female_snowd_old(styr,endyr,1,nlenm)
  matrix catch_female_rkd_new(styr,endyr,1,nlenm)
  matrix catch_female_rkd_old(styr,endyr,1,nlenm)
  matrix catch_trawl_female(styr,endyr,1,nlenm)
  matrix catch_trawl_male(styr,endyr,1,nlenm)
  
  matrix catch_male_ret(styr,endyr,1,nlenm)                         // Catch-in-numbers (males, retained)
  matrix catch_lmale(styr,endyr,1,nlenm)                            // Catch-in-numbers (males, total)
  3darray catch_discp(1,2,styr,endyr,1,nlenm)                       // Female discard catch (in mass)
  
  3darray natlength(1,2,styr,endyr,1,nlenm)                         // Total numbers by sex, length, and year
  3darray natlength_iold(1,2,styr,endyr,1,nlenm)                    // Immature old-shell numbers by sex, length, and year
  3darray natlength_inew(1,2,styr,endyr,1,nlenm)                    // Immature new-shell numbers by sex, length, and year
  3darray natlength_mold(1,2,styr,endyr,1,nlenm)                    // Mature old-shell numbers by sex, length, and year
  3darray natlength_mnew(1,2,styr,endyr,1,nlenm)                    // Mature new-shell numbers by sex, length, and year
  4darray natlength_mold_age(1,2,styr,endyr,1,nages,1,nlenm)        // Age- and length-structure
  3darray natlength_old(1,2,styr,endyr,1,nlenm)                     // Old-shell numbers by sex, length, and year
  3darray natlength_new(1,2,styr,endyr,1,nlenm)                     // New-shell numbers by sex, length, and year
  3darray natlength_i(1,2,styr,endyr,1,nlenm)                       // Immature numbers by sex, length, and year
  3darray natlength_mat(1,2,styr,endyr,1,nlenm)                     // Mature numbers by sex, length, and year
  3darray natl_new_fishtime(1,2,styr,endyr,1,nlenm)                 // Numbers-at-length (new shell)
  3darray natl_old_fishtime(1,2,styr,endyr,1,nlenm)                 // Numbers-at-length (old shell)
  3darray natl_inew_fishtime(1,2,styr,endyr,1,nlenm)                // Numbers-at-length (immature new shell)
  3darray natl_iold_fishtime(1,2,styr,endyr,1,nlenm)                // Numbers-at-length (immature new shell)
  3darray natl_mnew_fishtime(1,2,styr,endyr,1,nlenm)                // Numbers-at-length (mmature old shell)
  3darray natl_mold_fishtime(1,2,styr,endyr,1,nlenm)                // Numbers-at-length (mmature old shell)
  
  3darray len_len(1,2,1,nlenm,1,nlenm)                              // length to length growth array
  matrix moltp(1,2,1,nlenm)                                         // molting probabilities for female, male by length bin 
  matrix moltp_mat(1,2,1,nlenm)                                     // molting probs for mature female, male by length bin
  matrix mean_length(1,2,1,nlenm)                                   // Predicted post-moult sizes
  vector rec_len(1,nlenm)                                           // Recruitment length
  matrix rec_dev(1,2,styr,endyr-1)                                  // Recruitment deviations
  
  vector predpop_sexr(styr,endyr)                                   // Population sex-ratio - output
  
  vector fmort(styr,endyr)                                          // Directed (male) F
  vector fmortdf(styr,endyr)                                        // Directed (female) F
  vector fmortt(styr,endyr)
  vector fmortd1_snow(1992,endyr)
  vector fmortd1_rk(1,nobs_discardc-1)
  vector fmortd_snow(styr,endyr)
  vector fmortd_rk(styr,endyr)
  number breg_snow
  number breg_rk
                                                                    // Trawl F
  3darray F(1,2,styr,endyr,1,nlenm)                                 // Directed F on males
  3darray F_ret(1,2,styr,endyr,1,nlenm)                             // Retained F 
  matrix  Fdiscf(styr,endyr,1,nlenm)                                // Female target discards
  3darray Fdisct(1,2,styr,endyr,1,nlenm)
  3darray Fdisc_snow(1,2,styr,endyr,1,nlenm)
  3darray Fdisc_rk(1,2,styr,endyr,1,nlenm)
                                                                    // Trawl discards
  4darray Smat(1,2,1,2,styr,endyr,1,nlenm)                          // Survival 

  4darray effn_srv1(1,2,1,2,1,2,styr,endyr)                         // Survey effetive sample sizes
  matrix effn_fish_ret(1,2,styr,endyr)                              // Effective sample sizes:
  matrix effn_fish_tot(1,2,styr,endyr)                              // Effective sample sizes:
  
  // Offsets
  vector offset(1,8)
                                                  
  vector Fout(1,27);                                                // LIkelihoods and penalties
  
  // Penalties
  number nat_penalty
  number penal_rec                                                  // Recruitment
  number initsmo_penal                                              // Initial size-structure
  number initnum_penal                                              // Low size penalty
  number fpen                                                       // Penalties (misc)
  vector fmort_pen(1,7)                                             // vector of fmort dev penalties
  number sel_50m_penal                                              // Penalties on selectivity
  number af_penal                                                   // Prior on af 
  number srv3q_penalty
  number am_penal                                                   // Prior on am
  number bf_penal                                                   // Prior on bf
  number bm_penal                                                   // Prior on bm
  number penal_sexr                                                 // Penalty of sex ratio of recruitment            
  
  // Likelihood components
  number catch_like1                                                // catches - 1 
  number catch_like2                                                // catches - 2
  number catch_liket                                                // catches - 3
  number catch_likef                                                // catches - 4
  number catch_likes            // catches snow crab
  number catch_liker                                                // catches red king crab
  vector len_like(1,8)                                              // length data
  number cpue_like                                                  // cpue data
  number largemale_like                                             // Large males (why - AEP)
  number surv_like                                                  // Survey biomass data
  number like_initn
  number surv_like_nowt                                             // Surveys (output)
  3darray len_like_srv(1,2,1,2,1,2)                                 // Summary of survey likelihood
  
  vector emspbio_matetime(styr,endyr)                               // Spawning biomass stuff
  vector efspbio_matetime(styr,endyr)
  vector mspbio_matetime(styr,endyr)
  vector mspbio_old_matetime(styr,endyr)
  vector fspbio_new_matetime(styr,endyr)
  vector efspbio_new_matetime(styr,endyr)
  vector fspnum_new_matetime(styr,endyr)
  vector fspbio_matetime(styr,endyr)
  vector mspbio_fishtime(styr,endyr)
  vector fspbio_fishtime(styr,endyr)
  vector emspnum_old_matetime(styr,endyr)
  vector mspnum_matetime(styr,endyr)
  vector efspnum_matetime(styr,endyr)
  
  matrix maturity_est(1,2,1,nlenm)                                  // Maturity-at-length
  vector cpue_pred(styr,endyr)                                      // Predicted CPUE

  // Outputs (not in the likelihood function)
  vector num_males_gt101(styr,endyr)
  vector bio_males_gt101(styr,endyr)
  vector pred_catch_gt101(styr,endyr)
  vector pred_catch_no_gt101(styr,endyr)
  matrix pred_tmp(1,4,1,nobs_discardc)
  vector tmpp1(1,nlenm)
  vector tmpp2(1,nlenm)
  vector tmpp3(1,nlenm)
  vector tmpp4(1,nlenm)
  number like_mat
  
  sdreport_vector fspbios(1974,endyr)                                // Sd_report stuff
  sdreport_vector mspbios(1974,endyr)
  sdreport_vector legal_malesd(1974,endyr)
  sdreport_vector rec_early_sd(styr,1973)
  sdreport_vector recf_sd(1974,endyr-1)
  sdreport_vector recm_sd(1974,endyr-1)
  sdreport_number depletion
  objective_function_value f

 LOCAL_CALCS
   cout << "Phase: Moulting probabilities:       " << phase_moltingp << endl;
   cout << "Phase: Logistic selectivity pattern: " << phase_logistic_sel << endl;
   cout << "Phase: Dome-shaped selectivity:      " << phase_fishsel << endl;
   cout << "Phase: Survey selectivity #1         " << survq_phase << endl;
   cout << "Phase: Survey selectivity #2         " << survsel_phase << endl;
   cout << "Maturity switch:                     " << maturity_switch << endl;
   cout << "Growth switch:                       " << growth_switch << endl;
 END_CALCS

//  ========================================================================

PRELIMINARY_CALCS_SECTION
  int mat,shell,sex,l,ll,j,i,k,ii;  
  float ratio,sumtrawl,sumfishdiscf,sumPropn,sumfishdiscm,sumfishret,vall;

  // Specify initial values for the initial numbers
//  mnatlen_styr = log(10*(obs_p_srv1_lend(1,1,2,1)+obs_p_srv1_lend(2,1,2,1)+1e-02));
//  mnatlen_styr(2) = log(10*(obs_p_srv1_lend(1,2,2,1)+obs_p_srv1_lend(2,2,2,1)+1e-02));
//  for(j=1;j<=20;j++)
//   {
    //females first year of data is 1974 - fourth year of length data
//    fnatlen_styr(j) = log(10*(obs_p_srv1_lend(1,1,1,4,j)+obs_p_srv1_lend(2,1,1,4,j)+1e-02));
//    fnatlen_styr(2,j) = log(10*(obs_p_srv1_lend(1,2,1,4,j)+obs_p_srv1_lend(2,2,1,4,j)+1e-02));
//   }
// cout << "to here " << endl;
    
  // use logistic maturity curve for new shell males instead of fractions by year if switch>0
  // this would be for initial population not probability of moving to mature
  if(maturity_switch > 0)
  {
    maturity_average(2) = maturity_logistic;
  }
  log_input(catch_numbers);
//Compute offset for multinomial
  offset.initialize();
  for (i=1; i <= nobs_fish_discm; i++)
  {
    sumfishdiscm = sum(obs_p_fish_discmd(1,i)) + sum(obs_p_fish_discmd(2,i));
    for (j=1; j<=nlenm; j++)
    {
      obs_p_fish_tot(1,i,j)=obs_p_fish_discmd(1,i,j)/sumfishdiscm;
      obs_p_fish_tot(2,i,j)=obs_p_fish_discmd(2,i,j)/sumfishdiscm;
    }   
  }  
  // Offset for retained length-frequency
  for (i=1; i <= nobs_fish; i++)
  {
    sumfishret = sum(obs_p_fish_retd(1,i)) + sum(obs_p_fish_retd(2,i));
    obs_p_fish_ret(1,i) = obs_p_fish_retd(1,i)*fraction_new_error;
    obs_p_fish_ret(2,i) = obs_p_fish_retd(2,i)+(1.-fraction_new_error)*obs_p_fish_retd(1,i);
    for (j=1; j<=nlenm; j++)
    {
      for (shell=1; shell<=2; shell++)
      {
        obs_p_fish_ret(shell,i,j) = obs_p_fish_ret(shell,i,j)/sumfishret;
      }
      sumPropn = obs_p_fish_ret(1,i,j) + obs_p_fish_ret(2,i,j);
      offset(1) -= nsamples_fish(1,i)*sumPropn*log(p_const+sumPropn);
    }
  }  
  for (i=1; i <= nobs_fish_discm; i++)
  {
    for (j=1; j<=nlenm; j++)
    {
      sumPropn = obs_p_fish_tot(1,i,j)+obs_p_fish_tot(2,i,j);
      offset(2) -= nsamples_fish_discm(1,i)*sumPropn*log(p_const+sumPropn);
    } 
  }

  // Female discards 
  for (i=1; i <= nobs_fish_discf; i++)
  {
    sumfishdiscf = sum(obs_p_fish_discfd(i));
    for (j=1; j<=nlenm; j++)
    {
      obs_p_fish_discf(i,j) = obs_p_fish_discfd(i,j) / sumfishdiscf;
      offset(3) -= nsamples_fish_discf(i)*obs_p_fish_discf(i,j)*log(p_const+obs_p_fish_discf(i,j));
    }
  }  
  // Trawl discards
  for (i=1;i<=nobs_trawl;i++)
  {
    sumtrawl = 0;
    for(sex=1;sex<=2;sex++) 
      sumtrawl +=sum(obs_p_trawld(sex,i));
    for (j=1; j<=nlenm; j++)
    {
      obs_p_trawl(1,i,j) = obs_p_trawld(1,i,j) / sumtrawl;
      obs_p_trawl(2,i,j) = obs_p_trawld(2,i,j) / sumtrawl;
      offset(5)-=nsamples_trawl(1,i)*obs_p_trawl(1,i,j)*log(p_const +obs_p_trawl(1,i,j));
      offset(5)-=nsamples_trawl(2,i)*obs_p_trawl(2,i,j)*log(p_const+obs_p_trawl(2,i,j));
    }
  }
  //trawl catch biomass
  for (i=1;i<=nobs_trawl_c;i++){
    obs_catcht_biom(yrs_trawl_c(i))=catch_trawld(i);
  }
  // snow female discards
  for (i=1;i<=nobs_snowfish_discf;i++)
  {
    sumtrawl = 0.0;
    sumtrawl =sum(obs_p_snowfish_discf(i));
    for (j=1; j<=nlenm; j++)
    {
      obs_p_snow(1,i,j) = obs_p_snowfish_discf(i,j) / sumtrawl;
      offset(6)-=nsamples_snowfish_discf(i)*obs_p_snow(1,i,j)*log(p_const+obs_p_snow(1,i,j));
    }
  }

  // snow male discards  
  for (i=1;i<=nobs_snowfish_discm;i++)
  {
    sumtrawl = 0;
//  sex here is shell condition
    for(sex=1;sex<=2;sex++){ sumtrawl +=sum(obs_p_snowfish_discm(sex,i));}
// new and old shell together
    for (j=1; j<=nlenm; j++)
    {
      obs_p_snow(2,i,j) = (obs_p_snowfish_discm(1,i,j)+obs_p_snowfish_discm(2,i,j)) / sumtrawl;
      offset(6)-=nsamples_snowfish_discm(1,i)*obs_p_snow(2,i,j)*log(p_const+obs_p_snow(2,i,j));
    }
  }
  // red king female discards
  for (i=1;i<=nobs_rkfish_discf;i++)
  {
    sumtrawl = 0;
    sumtrawl +=sum(obs_p_rkfish_discf(i));
    for (j=1; j<=nlenm; j++)
    {
      obs_p_rk(1,i,j) = obs_p_rkfish_discf(i,j) / sumtrawl;
      offset(7)-= nsamples_rkfish_discf(i)*obs_p_rk(1,i,j)*log(p_const+obs_p_rk(1,i,j));
    }
  }
  // red king male discards
  for (i=1;i<=nobs_rkfish_discm;i++)
  {
    sumtrawl = 0;
    sumtrawl +=sum(obs_p_rkfish_discm(1,i)+obs_p_rkfish_discm(2,i));
    for (j=1; j<=nlenm; j++)
    {
      obs_p_rk(2,i,j) = (obs_p_rkfish_discm(1,i,j)+obs_p_rkfish_discm(2,i,j)) / sumtrawl;
      offset(8)-= nsamples_rkfish_discm(1,i)*obs_p_rk(2,i,j)*log(p_const+obs_p_rk(2,i,j));
    }
  }
  sumsrv.initialize();
  for(ll=1; ll<=nobs_srv1_length; ll++)
   for(mat=1;mat<=2;mat++)
    for(sex=1; sex<=2; sex++)
     for(j=1; j<=2; j++)
      sumsrv(ll) += sum(obs_p_srv1_lend(mat,sex,j,ll));

  // Survey data  
  for(mat=1; mat<=2; mat++) //maturity
   for(shell=1; shell<=2; shell++) //shell condition
    for(sex=1; sex<=2;sex++) //sex
     for (i=1; i <= nobs_srv1_length; i++)
      for (j=1; j<=nlenm; j++)
       {
        // only do new/old shell correction for mature crab
        if(mat < 2)
         obs_p_srv1_len1(mat,shell,sex,i,j) = obs_p_srv1_lend(mat,shell,sex,i,j)/sumsrv(i);
        else
         if( shell < 2)
          obs_p_srv1_len1(mat,shell,sex,i,j)=(obs_p_srv1_lend(mat,shell,sex,i,j)*fraction_new_error)/sumsrv(i);
         else
          obs_p_srv1_len1(mat,shell,sex,i,j)=(obs_p_srv1_lend(mat,shell,sex,i,j)+obs_p_srv1_lend(mat,1,sex,i,j)*(1.-fraction_new_error))/sumsrv(i);
       }

  // use logistic maturity curve for new and old shell male survey data if switch>0 instead of yearly samples
  // old shell already uses ok maturity curve (AEP only applies to OLD SHELL?)
  if(maturity_switch > 0)
   for(i=1; i <= nobs_srv1_length; i++)
    {
     tmps = (obs_p_srv1_len1(1,2,2,i)+obs_p_srv1_len1(2,2,2,i));
     obs_p_srv1_len1(2,2,2,i) = elem_prod(maturity_old_average(2),tmps);
     obs_p_srv1_len1(1,2,2,i) = elem_prod(1.0-maturity_old_average(2),tmps);
    }
   
  // Store results
  obs_p_srv1_len(1) = obs_p_srv1_len1(1);
  obs_p_srv1_len(2) = obs_p_srv1_len1(2);

  // for maturity and shell condition together in survey length comp fits
  offset(4) = 0.0;
  for(sex=1; sex<=2;sex++) //sex
   for (i=1; i <= nobs_srv1_length; i++)
    for (j=1; j<=nlenm; j++)
     { 
      vall = obs_p_srv1_len(1,1,sex,i,j)+obs_p_srv1_len(1,2,sex,i,j);
      offset(4) -= nsamples_srv1_length(2,1,sex,i)*vall*log(vall+p_const);
      vall = obs_p_srv1_len(2,1,sex,i,j)+obs_p_srv1_len(2,2,sex,i,j);
      offset(4) -= nsamples_srv1_length(2,2,sex,i)*vall*log(vall+p_const);
     } 

  cout << offset << endl;  
  log_input(offset);
                            
  obs_srv1_num.initialize();
  obs_srv1_biom.initialize();
  obs_srv1_bioms.initialize();
  obs_srv1_spbiom.initialize();
  obs_srv1_spnum.initialize();
  for(i=1;i<=nobs_srv1;i++) obs_srv1t(yrs_srv1(i)) = obs_srv1(i);
  
  // Compute survey biomass
  for(mat=1;mat<=2;mat++)  //maturity status
   for(l=1;l<=2;l++)  //shell condition
    for(sex=1;sex<=2;sex++)  //sex
     for (i=1; i <= nobs_srv1_length; i++)
      {
       obs_srv1_num(sex,yrs_srv1_length(i)) += obs_p_srv1_len(mat,l,sex,i)*obs_srv1t(yrs_srv1_length(i));
       if(sex<2)
        {
         obs_srv1_bioms(sex,yrs_srv1_length(i)) += obs_p_srv1_len(mat,l,sex,i)*obs_srv1t(yrs_srv1_length(i))*wtf(mat);
         obs_srv1_biom(yrs_srv1_length(i))    += obs_p_srv1_len(mat,l,sex,i)*obs_srv1t(yrs_srv1_length(i))*wtf(mat);
        }
       else
        {
         obs_srv1_bioms(sex,yrs_srv1_length(i)) += obs_p_srv1_len(mat,l,sex,i)*obs_srv1t(yrs_srv1_length(i))*wtm;
         obs_srv1_biom(yrs_srv1_length(i))    += obs_p_srv1_len(mat,l,sex,i)*obs_srv1t(yrs_srv1_length(i))*wtm;
        }        

       //  sum to get mature biomass by sex (AEP index is mature animals only?)
       if(mat>1)
        {
         if(sex<2)
          obs_srv1_spbiom(sex,yrs_srv1_length(i)) += obs_p_srv1_len(mat,l,sex,i)*obs_srv1t(yrs_srv1_length(i))*wtf(mat);
         else
          obs_srv1_spbiom(sex,yrs_srv1_length(i)) += obs_p_srv1_len(mat,l,sex,i)*obs_srv1t(yrs_srv1_length(i))*wtm;
         
         obs_srv1_spnum(l,sex,yrs_srv1_length(i)) += sum(obs_p_srv1_len(mat,l,sex,i)*obs_srv1t(yrs_srv1_length(i)));
        }
      }

  // Number of large males
  obs_lmales.initialize();
  obs_lmales_bio.initialize();
  for(i=1;i<=nobs_srv1_length;i++)
   {
    // take 1/2 of the 100-104 bin, 
    obs_lmales(i) = 0.5*obs_srv1_num(2,yrs_srv1_length(i),23);
    obs_lmales_bio(i) = obs_lmales(i)*wtm(23);
    for(j=24;j<=nlenm;j++)
     {
      obs_lmales(i) += obs_srv1_num(2,yrs_srv1_length(i),j);
      obs_lmales_bio(i) += obs_srv1_num(2,yrs_srv1_length(i),j)*wtm(j);
     }
   }
  // average catch (maes and females) by trawl
 // avgpf=0; avgpm=0;
 // for(i=1; i <= nobs_trawl; i++)
 //  { avgpf+=obs_p_trawl(1,i); avgpm+=obs_p_trawl(2,i); }
//  avgpf=avgpf/nobs_trawl;
//  avgpm=avgpm/nobs_trawl;
 
  // Total trawl catch in mass
//  for(i=styr;i<=(yrs_trawl(1)-1);i++)
//   obs_catcht_biom(i) = avgpf*catch_trawl(i)*wtf(2) + avgpm*catch_trawl(i)*wtm;
//  cout<<"to beg disc biomass "<<endl;
  // Multiply total numbers by avearge weight
  obs_catchdf_biom.initialize();
   avgp.initialize();
  for(i=1;i<=nobs_fish_catchf;i++)
  {
    obs_catchdm_biom(yrs_fish_catchf(i)) = (catch_odisc(2,i));
    obs_catchdf_biom(yrs_fish_catchf(i)) = catch_odisc(1,i);
    obs_catchtot_biom(yrs_fish_catchf(i)) = obs_catchdm_biom(yrs_fish_catchf(i))+catch_ret(yrs_fish_catchf(i));}
  
  // Historical approximation to observed catch biomass
  // Compute the moulting probabilities
  get_moltingp();
  // estimate growth function
  get_growth();
  // Set maturity
  get_maturity();
  cout<<"end prelim calcs"<<endl;
// ============================================================================
PROCEDURE_SECTION
  // Update growth (if the parameters are being estimated)
  if (active(moltp_af) || active(moltp_bf) || active(moltp_am) || active(moltp_bm) || active(moltp_ammat) || active(moltp_bmmat)) get_moltingp();
  // growth estimated in prelimn calcs if growth parameters estimated in the model
  // then will redo growth matrix, otherwise not
  if(active(am1) || active(bm1) || active(af1) || active(bf1) || active(growth_beta)) get_growth();
   get_maturity();
   get_selectivity();
   get_mortality();
   get_numbers_at_len();
   get_catch_at_len();
   evaluate_the_objective_function();
// ----------------------------------------------------------------------
FUNCTION WriteMCMC
 post<<
// srv1_slope <<","<<
// srv1_sel50 <<","<<
 fish_slope_mn <<","<<
 fish_sel50_mn <<","<<
// fish_fit_slope_mn <<","<<
// fish_fit_sel50_mn <<","<<
 fish_disc_slope_f <<","<<
 fish_disc_sel50_f <<","<<
 //fish_disc_slope_tf <<","<<
 //fish_disc_sel50_tf <<","<<
 endl;

// --------------------------------------------------------------------------

FUNCTION get_maturity
  int j;
  if(active(matestm))
  {
    maturity_est(1)(1,16) = mfexp(matestf);
    maturity_est(1)(17,nlenm) = 1.;
    maturity_est(2) = mfexp(matestm);
  }
  else
  {    
    maturity_est(1) = maturity_average(1);
    maturity_est(2) = maturity_average(2);
  }
// --------------------------------------------------------------------------

FUNCTION get_growth
  int ilen,il2,sex;
  dvariable devia, recsum, alpha_rec,alpha,lensum;
  dvariable growinc_67, growinc_90;
  len_len.initialize();
//  growinc_90 = am2 + bm2 * 90.0 - 90.0;
//  growinc_67 = am1 + bm1 * 67.5 - 67.5;
  for(ilen=1;ilen<=nlenm;ilen++)
  { 
    //linear growth curve
    if(growth_switch==1)
    {
      mean_length(1,ilen)= mfexp(af1)* pow(length_bins(ilen),bf1);
     if(ilen<17){    
        mean_length(2,ilen)= mfexp(am1)* pow(length_bins(ilen),bm1);
      }
      else{
        mean_length(2,ilen)= mfexp(am1)* pow(length_bins(ilen),bm1);
      }
    }
    // growth_switch==2 constant after 90mm
//    if(growth_switch==2)
//     {  
//      if(ilen<14)
//       {
//        mean_length(1,ilen)= af + bf * length_bins(ilen);
//        mean_length(2,ilen)= am1 + bm1 * length_bins(ilen);
//       }
//     else
//      {
//       mean_length(1,ilen)= mean_length(1,13)-length_bins(13)+length_bins(ilen);
//       mean_length(2,ilen)= growinc_90 +length_bins(ilen);
//      }
//     }
    // this is growth from Otto(1995) to 67 mm then constant after
//    if(growth_switch==3)
//     {
//      mean_length(1,ilen)= af + bf * length_bins(ilen); 
//      if(ilen < 10)
//       mean_length(2,ilen)= am1 + bm1 * length_bins(ilen);
 //     else
 //      mean_length(2,ilen)= growinc_67 +length_bins(ilen);
//     }
   }
// using Gamma function for transition matrix
// devia is the bounds of growth bins to evaluate
// the gamma function (x) in prop = integral(i1 to i2) g(x|alpha,beta) dx
// alpha and growth_beta are parameters 
// alpha is the mean growth increment per molt for some premolt length class
// alpha = mean growth increment per molt divided by beta
// beta is the shape parameter - larger beta - more variance 
 
  for (sex=1;sex<=2;sex++)
  {
    for(ilen=1;ilen<=nlenm;ilen++)
    {
     // subract the 2.5 from the midpoint of the length bin to get the lower bound
      alpha = (mean_length(sex,ilen)-(length_bins(ilen)-2.5))/growth_beta(sex);
      lensum = 0;
// this statement non-truncated growth transition
// for(il2=ilen;il2<=nlenm;il2++)
//    truncate growth transition to max=10 bins
      for(il2=ilen;il2<=ilen+min(10,nlenm-ilen);il2++)
      {
          devia = length_bins(il2)+2.5-length_bins(ilen);
          len_len(sex,ilen,il2) = pow(devia,(alpha-1.))*exp(-devia/growth_beta(sex));
//        len_len(sex,ilen,il2) = mfexp( log(devia)*(alpha-1.) - devia/growth_beta(sex));
          lensum += len_len(sex,ilen,il2);
      }  
      len_len(sex,ilen) /= sum(len_len(sex,ilen));
    }
  }
 // Fraction recruiting
  recsum=0.0;
  alpha_rec=alpha1_rec/beta_rec;
  for(ilen=1;ilen<=nlenm;ilen++)
  {
    devia = length_bins(ilen)+2.5-length_bins(1);
//    rec_len(ilen) =  mfexp( log(devia)*(alpha_rec-1.) - devia/beta_rec);
    rec_len(ilen) =  pow(devia,alpha_rec-1.)*exp(-devia/beta_rec);
    recsum += rec_len(ilen);
  }
  //standardize so each row sums to 1.0
  for(ilen=1;ilen<=nlenm;ilen++)
    rec_len(ilen) = rec_len(ilen)/recsum;
// -------------------------------------------------------------------------
FUNCTION get_moltingp
  int j;
  //assuming a declining logistic function
  for(j=1;j<=nlenm;j++)
  {
    // logistic molting females then males
    moltp(1,j)=1-(1./(1.+mfexp(-1.*moltp_af*(length_bins(j)-moltp_bf))));
    moltp(2,j)=1-(1./(1.+mfexp(-1.*moltp_am*(length_bins(j)-moltp_bm))));

    // set molting prob for mature females at 0.0
    moltp_mat(1,j)=0.0;
    // molting probability can be one (or estimated)
    if(phase_moltingp > 0)
     moltp_mat(2,j) = 1-(1./(1.+mfexp(-1.*moltp_ammat*(length_bins(j)-moltp_bmmat))));
    else
     moltp_mat(2,j) = 0.0;
  }

// -------------------------------------------------------------------------

FUNCTION get_selectivity
  int iy,j,ii;
  dvariable maxsel_fish;
  dvariable tmp2;
  ii=1;
  // logistic selectivity curves
  if(active(log_sel50_dev_mn))
  {
    for(iy=styr;iy<=endyr;iy++)
    {
     // length at 50% selectivity
      if(iy < 1980 || (iy > 1984 && iy < 1988) || (iy > 1996 && iy < 2005  || iy > 2009))
        fish_sel50_mn(iy)=mfexp(log_avg_sel50_mn);
      else
      {
        fish_sel50_mn(iy)=mfexp(log_avg_sel50_mn+log_sel50_dev_mn(ii));
        ii=ii+1;
      }
      //logistic selectivity curve
      for (j=1;j<=nlenm;j++)
      { 
        sel(1,iy,j)=1./(1.+mfexp(-1.*fish_slope_mn*(length_bins(j)-fish_sel50_mn(iy))));
        if(iy < 1993)    
          sel_ret(1,iy,j) = 1./(1.+mfexp(-1.*fish_fit_slope_mn1*(length_bins(j)-fish_fit_sel50_mn1)));
        else
          sel_ret(1,iy,j) = 1./(1.+mfexp(-1.*fish_fit_slope_mn2*(length_bins(j)-fish_fit_sel50_mn2)));
        //for dome shaped add this part
        if(phase_fishsel > 0)
        {
          tmp2=1./(1.+mfexp(fish_slope_mn2*(length_bins(j)-fish_sel50_mn2)));
          sel(1,iy,j)=sel(1,iy,j)*tmp2;
        }
 
        // set new and old sel same
        sel(2,iy,j)=sel(1,iy,j);
        sel_ret(2,iy,j)=sel_ret(1,iy,j);
      }
    }
  }
 
  for(iy=styr;iy<=endyr;iy++)
  {
     if (iy<1991)
       sel(1,iy)=1./(1.+mfexp(-fish_slope_1*(length_bins-mean(mfexp(log_avg_sel50_3+log_sel50_dev_3(1,6))))));
     if(iy>1990 && iy < 1997)
       sel(1,iy)=1./(1.+mfexp(-fish_slope_1*(length_bins-mfexp(log_avg_sel50_3+log_sel50_dev_3(iy-1990)))));
     if(iy>1996 && iy<2005)
       sel(1,iy)=1./(1.+mfexp(-fish_slope_yr_3*(length_bins-mfexp(log_avg_sel50_3))));
     
     if(iy>2004 && iy<(endyr-1))
       sel(1,iy)=1./(1.+mfexp(-fish_slope_yr_3*(length_bins-mfexp(log_avg_sel50_3+log_sel50_dev_3(iy-1998)))));
     if(iy>(endyr-2))
       sel(1,iy)=1./(1.+mfexp(-fish_slope_yr_3*(length_bins-mfexp(log_avg_sel50_3))));
     if(iy < 1991)    
       sel_ret(1,iy) = 1./(1.+mfexp(-1.*fish_fit_slope_mn1*(length_bins-fish_fit_sel50_mn1)));
     else
       sel_ret(1,iy) = 1./(1.+mfexp(-1.*fish_fit_slope_mn2*(length_bins-fish_fit_sel50_mn2)));
    //for dome shaped add this part
    if(phase_fishsel > 0)
    {
      sel(1,iy)= elem_div( sel(1,iy), (1.+mfexp(fish_slope_mn2*(length_bins-fish_sel50_mn2))) );
    }
    // set new and old sel same
    sel(2,iy)     = sel(1,iy);
    sel_ret(2,iy) = sel_ret(1,iy);
  }
// female discards ascending logistic curve 
    sel_discf=1./(1.+mfexp(-1.*fish_disc_slope_f*(length_bins-fish_disc_sel50_f)));
//      sel_disc_snow(1,j)=1./(1.+mfexp(-1.*snowfish_disc_slope_f*(length_bins-snowfish_disc_sel50_f)));
      
//  snow fishery selectivity for 3 time periods, #1 (1989-1996), #2 (1997-2004) and #3 (2005-P)      
    sel_disc_snow(1,1)=1./(1.+mfexp(-1.*snowfish_disc_slope_f_1*(length_bins-snowfish_disc_sel50_f_1))); 
    sel_disc_snow(2,1)=1./(1.+mfexp(-1.*snowfish_disc_slope_f_2*(length_bins-snowfish_disc_sel50_f_2))); 
    sel_disc_snow(3,1)=1./(1.+mfexp(-1.*snowfish_disc_slope_f_3*(length_bins-snowfish_disc_sel50_f_3))); 
     
//snow males dome shaped
//      sel_disc_snow(2)=(1./(1.+mfexp(-1.*snowfish_disc_slope_m*(length_bins-snowfish_disc_sel50_m)))) *(1./(1.+mfexp(snowfish_disc_slope_m2*(length_bins-snowfish_disc_sel50_m2))));

//  snow fishery selectivity for 3 time periods, #1 (1989-1996), #2 (1997-2004) and #3 (2005-P)      
    sel_disc_snow(1,2)= elem_div(1./(1.+mfexp(-1.*snowfish_disc_slope_m_1*(length_bins-snowfish_disc_sel50_m_1))), (1.+mfexp(snowfish_disc_slope_m2_1*(length_bins-snowfish_disc_sel50_m2_1))));
    sel_disc_snow(2,2)= elem_div(1./(1.+mfexp(-1.*snowfish_disc_slope_m_2*(length_bins-snowfish_disc_sel50_m_2))), (1.+mfexp(snowfish_disc_slope_m2_2*(length_bins-snowfish_disc_sel50_m2_2))));
    sel_disc_snow(3,2)= elem_div(1./(1.+mfexp(-1.*snowfish_disc_slope_m_3*(length_bins-snowfish_disc_sel50_m_3))), (1.+mfexp(snowfish_disc_slope_m2_3*(length_bins-snowfish_disc_sel50_m2_3))));
      
//  red king crab fishery selectivity for 3 time periods, #1 (1989-1996), #2 (1997-2004) and #3 (2005-P)      
//  red fishery selectivity for 3 time periods, #1 (1989-1996), #2 (1997-2004) and #3 (2005-P)      
    sel_disc_rkc(1,1)=1./(1.+mfexp(-1.*rkfish_disc_slope_f1*(length_bins-rkfish_disc_sel50_f1))); 
    sel_disc_rkc(2,1)=1./(1.+mfexp(-1.*rkfish_disc_slope_f2*(length_bins-rkfish_disc_sel50_f2))); 
    sel_disc_rkc(3,1)=1./(1.+mfexp(-1.*rkfish_disc_slope_f3*(length_bins-rkfish_disc_sel50_f3))); 
       
 //  red fishery selectivity for 3 time periods, #1 (1989-1996), #2 (1997-2004) and #3 (2005-P)      
    sel_disc_rkc(1,2)=1./(1.+mfexp(-1.*rkfish_disc_slope_m1*(length_bins-rkfish_disc_sel50_m1))); 
    sel_disc_rkc(2,2)=1./(1.+mfexp(-1.*rkfish_disc_slope_m2*(length_bins-rkfish_disc_sel50_m2))); 
    sel_disc_rkc(3,2)=1./(1.+mfexp(-1.*rkfish_disc_slope_m3*(length_bins-rkfish_disc_sel50_m3))); 
       
 // trawl fishery selectivity
 
 //  trewl fishery selectivity for 3 time periods, #1 (1973-1987), #2 (1988-1996) and #3 (1997-P)
    sel_trawl(1,1)=1./(1.+mfexp(-1.*fish_disc_slope_tf1*(length_bins-fish_disc_sel50_tf1)));
    sel_trawl(1,2)=1./(1.+mfexp(-1.*fish_disc_slope_tm1*(length_bins-fish_disc_sel50_tm1)));

    sel_trawl(2,1)=1./(1.+mfexp(-1.*fish_disc_slope_tf2*(length_bins-fish_disc_sel50_tf2)));
    sel_trawl(2,2)=1./(1.+mfexp(-1.*fish_disc_slope_tm2*(length_bins-fish_disc_sel50_tm2)));
     
    sel_trawl(3,1)=1./(1.+mfexp(-1.*fish_disc_slope_tf3*(length_bins-fish_disc_sel50_tf3)));
    sel_trawl(3,2)=1./(1.+mfexp(-1.*fish_disc_slope_tm3*(length_bins-fish_disc_sel50_tm3)));

  for(iy=styr;iy<=endyr;iy++)
  {
    // maxsel_fish=max(sel(1,iy));
    maxsel_fish   = sel(1,iy,15);
    if(maxsel_fish<max(sel(2,iy))) 
      maxsel_fish = sel(2,iy,15);
      // maxsel_fish=max(sel(2,iy));
    sel(1,iy)     = sel(1,iy)/maxsel_fish;
    sel(2,iy)     = sel(2,iy)/maxsel_fish;
    sel_fit(1,iy) = elem_prod(sel_ret(1,iy),sel(1,iy));
    sel_fit(2,iy) = elem_prod(sel_ret(2,iy),sel(2,iy));
  }
  get_srv_selectivity();

FUNCTION get_srv_selectivity
    // somerton and otto curve for survey selectivities
    // need to get time-periods for selectivity
    if (survsel_som_phase<0)
      // This would need fixing...to work with tc version 2
      sel_srv_m(styr) = 1./(1.+sel_som(2)*mfexp(-1.*sel_som(3)*length_bins));
    else
    {
      for (int iyr=styr;iyr<yrs_q(1);iyr++)
      {
        srv_q_f(iyr) = srv_q_f_par;
        srv_q_m(iyr) = srv_q_m_par;
        sel_srv_f(iyr) = srv_q_f(iyr) /(1.+mfexp(-1.*log(19.)*(length_bins-srv_sel50_f)/(srv_sel95_f-srv_sel50_f)));
        sel_srv_m(iyr) = srv_q_m(iyr) /(1.+mfexp(-1.*log(19.)*(length_bins-srv_sel50_m)/(srv_sel95_m-srv_sel50_m)));
      }
      dvariable srvq_f_tmp ;
      dvariable srvq_m_tmp ;
      dvariable srv50_f_tmp ;
      dvariable srv50_m_tmp ;
      dvariable srv95_f_tmp ;
      dvariable srv95_m_tmp ;
      for (int iq=1;iq<=n_qs;iq++)
      {
        if (iq>1)
        {
          srvq_f_tmp = srv_q_f_par*mfexp(srv_q_f_dev(iq-1));
          srvq_m_tmp = srv_q_m_par*mfexp(srv_q_m_dev(iq-1));
          srv50_f_tmp = srv_sel50_f*mfexp(srv_sel50_f_dev(iq-1));
          srv50_m_tmp = srv_sel50_m*mfexp(srv_sel50_m_dev(iq-1));
          srv95_f_tmp = srv_sel95_f*mfexp(srv_sel95_f_dev(iq-1));
          srv95_m_tmp = srv_sel95_m*mfexp(srv_sel95_m_dev(iq-1));
        }
        else
        {
          srvq_f_tmp = srv_q_f_par;
          srvq_m_tmp = srv_q_m_par;
          srv50_f_tmp = srv_sel50_f;
          srv50_m_tmp = srv_sel50_m;
          srv95_f_tmp = srv_sel95_f;
          srv95_m_tmp = srv_sel95_m;
        }
        int period_end;
        if (iq<n_qs)
          period_end=yrs_q(iq+1);
        else
          period_end=endyr+1;
        // In between years of q's
        for (int iyr=yrs_q(iq);iyr<period_end;iyr++)
        {
          // Set q, then selectivity 
          srv_q_f(iyr) = srvq_f_tmp;
          srv_q_m(iyr) = srvq_m_tmp;
          sel_srv_f(iyr) = srv_q_f(iyr) / (1.+mfexp(-1.*log(19.)*(
                           length_bins-srv50_f_tmp)/(srv95_f_tmp-srv50_f_tmp)));
          sel_srv_m(iyr) = srv_q_m(iyr) / (1.+mfexp(-1.*log(19.)*(
                           length_bins-srv50_m_tmp)/(srv95_m_tmp-srv50_m_tmp)));
        }
      }
    }
// -------------------------------------------------------------------------

FUNCTION get_mortality
  int i,shell;
  int inc,ii;
  inc=0;
  M(1)= M_in(1)*Mmult_imat;
  M(2)= M_in(2)*Mmult_imat;
  M_matn(2)= M_matn_in(2)*Mmultm;
  M_mato(2)= M_mato_in(2)*Mmultm;
  M_matn(1)= M_matn_in(1)*Mmultf;
  M_mato(1)= M_mato_in(1)*Mmultf;
  
  //first year retained catch 1965(1966 fishery) no fishery 86, 87 or 1997 to 2004
  for(i =styr;i<=endyr;i++){
    if(i>1965){ 
      if(i==1985 || i==1986 || (i>1996 && i<2005) || i==2010 || i==2011){
        fmort(i)=0.0;
      }
      else{
        ii=1966+inc;
        fmort(i) = mfexp(log_avg_fmort+fmort_dev(ii));
        inc=inc+1;
      }
    }
    else{
      fmort(i)= 0.05;  //mean(fmort(1965,1979));
    }
  }
// fmortdf=mfexp(log_avg_fmortdf+fmortdf_dev); using overall F for females as well as males in directed fishery
  fmortt=0.0000001;
 //Fs in snow and rkc fishery are scalars need to multiply in projections by retained snow crab/average snow catch * F to get F.
  fmortd1_snow = mfexp(log_avg_fmortd_snow+fmortd_snow_dev);
// breg_snow=(sum(elem_prod(fmortd1_snow,snowcatch(1992,endyr)))-(sum(snowcatch(1992,endyr))*sum(fmortd1_snow)/(endyr-1992)))/(norm2(snowcatch(1992,endyr))-((sum(snowcatch(1992,endyr))*sum(snowcatch(1992,endyr)))/(endyr-1992)));
  breg_snow = mean(fmortd1_snow)/(mean(snowcatch(1992,endyr)));
  for(i=styr;i<=endyr;i++){
    if(i<1992)
    {
      if(i<1974)
        fmortd_snow(i)= 0.01;  //fmortd_snow(i)= breg_snow * sum(snowcatch(1974,1978))/5.0;
      else 
        fmortd_snow(i) = 0.0015*snowcatch(i);       //fmortd_snow(i) = breg_snow*snowcatch(i);
    }
    else
      fmortd_snow(i) = mfexp(log_avg_fmortd_snow+fmortd_snow_dev(i));
  }
  for(i=styr;i<=endyr;i++)
  {
    if(i>1972)
      fmortt(i) = mfexp(log_avg_fmortt+fmortt_dev(i));
    else
      fmortt(i) = mfexp(log_avg_fmortt);
      // fmortt(i) = mean(mfexp(log_avg_fmortt+fmortt_dev));
  }
// need to have the devs 1992 to present
  fmortd1_rk = mfexp(log_avg_fmortd_rk+fmortd_rk_dev);
// breg_rk=(sum(elem_prod(fmortd1_rk,rkccatch(1992,endyr)))-(sum(rkccatch(1992,endyr))*sum(fmortd1_rk)/(endyr-1992)))/(norm2(rkccatch(1992,endyr))-((sum(rkccatch(1992,endyr))*sum(rkccatch(1992,endyr)))/(endyr-1992)));
  breg_rk = mean(1.-exp(-fmortd1_rk))/(mean(rkccatch(yrs_discardc)));
  ii=1;
  for(i=styr;i<=endyr;i++)
  {
    if(i<1969){
      fmortd_rk(i)= 0.02; //breg_rk*mean(rkccatch(1969,1973));
    }
    else{
      if(i>1968 && i<1992)
      { 
        fmortd_rk(i) = -log(1.-breg_rk*rkccatch(i));
        // This is incorrect--not differentiable...
        // if(fmortd_rk(i)< 0.01) fmortd_rk(i)= 0.01;
        if(i>1983 && i <1986) fmortd_rk(i)=0.0;
      }
      else
      { 
        fmortd_rk(i) = mfexp(log_avg_fmortd_rk+fmortd_rk_dev(ii));
        if(i>1993 && i<1996)
        { 
          fmortd_rk(i)=0.0;
        }
        else{
          ii=ii+1;}
      }
    }
  }
  for (i=styr;i<=endyr;i++)
  {
//need to set F in directed fishery to 0.0 when was closed 1985-1986 and 1997-2004 and?
//set F in red king to 0 when closed 84-85 and 94-95
   //  have discard mort for females and males, fishing F for males only
   //   Fdisct(1,i)= sel_trawl(1)*fmortt(i);   
   //   Fdisct(2,i)= sel_trawl(2)*fmortt(i);
    Fdiscf(i)= sel_discf*fmort(i);
 
 // test on year for 3 trawl selectivity periods
    if (i<1987)
    {
      sel_trawl_use(1)=sel_trawl(1,1);
      sel_trawl_use(2)=sel_trawl(1,2);
    }
    if (i>1986 && i<1997)
    {
      sel_trawl_use(1)=sel_trawl(2,1);
      sel_trawl_use(2)=sel_trawl(2,2);
    }
    if (i>1996)
    {
      sel_trawl_use(1)=sel_trawl(3,1);
      sel_trawl_use(2)=sel_trawl(3,2);
    }
 
 // test on year for 3 snow selectivity periods
    if (i<1997)
    {
      sel_disc_snow_use(1)=sel_disc_snow(1,1);
      sel_disc_snow_use(2)=sel_disc_snow(1,2);
    }
    if (i>1996 && i<2005)
    {
      sel_disc_snow_use(1)=sel_disc_snow(2,1);
      sel_disc_snow_use(2)=sel_disc_snow(2,2);
    }
    if (i>2004)
    {
      sel_disc_snow_use(1)=sel_disc_snow(3,1);
      sel_disc_snow_use(2)=sel_disc_snow(3,2);
    }
 
    // test on year for 3 red selectivity periods
    if (i<1997)
    {
      sel_disc_rkc_use(1)=sel_disc_rkc(1,1);
      sel_disc_rkc_use(2)=sel_disc_rkc(1,2);
    }
    if (i>1996 && i<2005)
    {
      sel_disc_rkc_use(1)=sel_disc_rkc(2,1);
      sel_disc_rkc_use(2)=sel_disc_rkc(2,2);
    }
    if (i>2004)
    {
      sel_disc_rkc_use(1)=sel_disc_rkc(3,1);
      sel_disc_rkc_use(2)=sel_disc_rkc(3,2);
    }
      
    Fdisct(1,i)= sel_trawl_use(1)*fmortt(i);   
    Fdisct(2,i)= sel_trawl_use(2)*fmortt(i);   
    Fdisc_snow(1,i)= sel_disc_snow_use(1)*fmortd_snow(i);   
    Fdisc_snow(2,i)= sel_disc_snow_use(2)*fmortd_snow(i);   
    Fdisc_rk(1,i)= sel_disc_rkc_use(1)*fmortd_rk(i);   
    Fdisc_rk(2,i)= sel_disc_rkc_use(2)*fmortd_rk(i);   
 
    for(shell=1;shell<=2;shell++) //over new (shell=1) and old (shell=2) shell...
    { 
      F(shell,i) = sel(shell,i)*fmort(i);       
      F_ret(shell,i)=sel_fit(shell,i)*fmort(i);
      Smat(1,shell,i)=mfexp(-1.0*(Fdiscf(i) + Fdisct(1,i)+Fdisc_snow(1,i)+Fdisc_rk(1,i)));
      Smat(2,shell,i)=mfexp(-1.0*(F(shell,i)+ Fdisct(2,i)+Fdisc_snow(2,i)+Fdisc_rk(2,i)));
    } 
  }

// -------------------------------------------------------------------------

FUNCTION get_numbers_at_len
  int itmp, jk, sex, j, i;
  dvar_matrix tmpo(1,2,styr,endyr);
  dvariable tmpi,Surv1,Surv2,Surv3,Surv4,Surv5,Surv6;
  dvariable totSrvNum;
   
  natlength.initialize();
  natlength_inew.initialize();
  natlength_iold.initialize();
  natlength_mnew.initialize();
  natlength_mold.initialize();
  natlength_new.initialize();
  natlength_old.initialize();
  natlength_i.initialize();
  natlength_mat.initialize();
  mean_log_rec(1) = mean_log_rec1;
  mean_log_rec(2) = mean_log_rec1;
  rec_dev(1)(1974,endyr-1) = rec_devf;
  rec_dev(2)(1974,endyr-1) = rec_devf;
  rec_dev(1)(styr,1973) = rec_devf_early;
  rec_dev(2)(styr,1973) = rec_devf_early;

  //numbers at length from styr to endyr
  for(sex=1;sex<=2;sex++)
  {
    for (i=styr;i< endyr;i++)
    {
      Surv1 = mfexp(-(1-catch_midpt(i))*M(sex));
      Surv1 = mfexp(-M(sex));
      Surv2 = mfexp(-catch_midpt(i)*M(sex));
      Surv3 = mfexp(-(1-catch_midpt(i))*M_matn(sex));
      Surv3 = mfexp(-M_matn(sex));
      Surv4 = mfexp(-catch_midpt(i)*M_matn(sex));
      Surv5 = mfexp(-M_mato(sex));
      Surv6 = mfexp(-catch_midpt(i)*M_mato(sex));

      if(i==1983 && killem_switch==1)
      {
        Surv3 = mfexp(-(1-catch_midpt(i))*M_matn(sex)*mat_big(sex));
        Surv3 = mfexp(-M_matn(sex)*mat_big(sex));
        Surv4 = mfexp(-catch_midpt(i)*M_matn(sex)*mat_big(sex));
        Surv5 = mfexp(-M_mato(sex)*mat_big(sex));
        Surv6 = mfexp(-catch_midpt(i)*M_mato(sex)*mat_big(sex));
      }       
      if(i==styr)
      { 
        // add in recruits for next year
        // put all recruits in new shell immature
        natlength_inew(sex,i) += mfexp(mean_log_rec1_early+rec_devf_early(i))*rec_len*proprecn;
      }
      // Numbers advancing to new shell...
      dvar_vector tmp = Surv1*elem_prod(moltp(sex),elem_prod(Smat(sex,1,i),natlength_inew(sex,i)));
      natlength_new(sex,i+1) =  tmp * len_len(sex);
      
      dvar_vector tmpo = Surv1*elem_prod(moltp(sex),elem_prod(Smat(sex,2,i),natlength_iold(sex,i)));
      natlength_new(sex,i+1) +=  tmpo * len_len(sex);

      natlength_iold(sex,i+1) = Surv1*(elem_prod(Smat(sex,1,i),natlength_inew(sex,i)) + 
                                       elem_prod(Smat(sex,2,i),natlength_iold(sex,i))) - tmp-tmpo;
 
      dvar_vector tmpm = Surv3*elem_prod(moltp_mat(sex),elem_prod(Smat(sex,1,i),natlength_mnew(sex,i)));
      natlength_mnew(sex,i+1) =   tmpm * len_len(sex);

      dvar_vector tmpmo = Surv5*elem_prod(moltp_mat(sex),elem_prod(Smat(sex,2,i),natlength_mold(sex,i)));
      natlength_mnew(sex,i+1) +=   tmpmo * len_len(sex);
 
      natlength_mold(sex,i+1) = Surv3 * elem_prod(Smat(sex,1,i),natlength_mnew(sex,i)) + 
                                Surv5 * elem_prod(Smat(sex,2,i),natlength_mold(sex,i)) - tmpm-tmpmo;
      
      // this is for estimating the fraction of new shell that move to old shell to fit
      // the survey data that is split by immature and mature
      natlength_mnew(sex,i+1) += elem_prod(    maturity_est(sex),natlength_new(sex,i+1));
      natlength_inew(sex,i+1)  = elem_prod(1.0-maturity_est(sex),natlength_new(sex,i+1));
      // add in recruits for next year
      // put all recruits in new shell immature
      if(i<1974)          
        natlength_inew(sex,i+1) += mfexp(mean_log_rec1_early+rec_devf_early(i))*rec_len*proprecn;
      else
        natlength_inew(sex,i+1) += mfexp(mean_log_rec(sex)+rec_dev(sex,i))*rec_len*proprecn;

      natlength_new(sex,i+1)     = natlength_inew(sex,i+1) + natlength_mnew(sex,i+1);
      natlength_old(sex,i+1)     = natlength_mold(sex,i+1) + natlength_iold(sex,i+1);
      natlength_mat(sex,i+1)     = natlength_mnew(sex,i+1) + natlength_mold(sex,i+1);
      natlength_i(sex,i+1)       = natlength_inew(sex,i+1) + natlength_iold(sex,i+1);
      natlength(sex,i+1)         = natlength_mat(sex,i+1)  + natlength_i(sex,i+1);
    }
  }
  for (i=styr;i<=endyr;i++)
  {
    for (sex=1;sex<=2;sex++)
    {
      if(i==1983 && killem_switch==1)
      {
        natl_inew_fishtime(sex,i) = mfexp(-catch_midpt(i)*M(sex))*natlength_inew(sex,i);
        natl_iold_fishtime(sex,i) = mfexp(-catch_midpt(i)*M(sex))*natlength_iold(sex,i);
        natl_mnew_fishtime(sex,i) = mfexp(-catch_midpt(i)*M_matn(sex)*mat_big(sex))*natlength_mnew(sex,i);
        natl_mold_fishtime(sex,i) = mfexp(-catch_midpt(i)*M_mato(sex)*mat_big(sex))*natlength_mold(sex,i);
        natl_new_fishtime(sex,i)  = natl_inew_fishtime(sex,i)+natl_mnew_fishtime(sex,i); 
        natl_old_fishtime(sex,i)  = natl_iold_fishtime(sex,i)+natl_mold_fishtime(sex,i);
      }
      else
      {
        natl_inew_fishtime(sex,i) = mfexp(-catch_midpt(i)*M(sex))*natlength_inew(sex,i);
        natl_iold_fishtime(sex,i) = mfexp(-catch_midpt(i)*M(sex))*natlength_iold(sex,i);
        natl_mnew_fishtime(sex,i) = mfexp(-catch_midpt(i)*M_matn(sex))*natlength_mnew(sex,i);
        natl_mold_fishtime(sex,i) = mfexp(-catch_midpt(i)*M_mato(sex))*natlength_mold(sex,i);
        natl_new_fishtime(sex,i)  = natl_inew_fishtime(sex,i)+natl_mnew_fishtime(sex,i); 
        natl_old_fishtime(sex,i)  = natl_iold_fishtime(sex,i)+natl_mold_fishtime(sex,i);
      }
    }
    totn_srv1(1,i) = natlength(1,i)*srv_q_f(i)*sel_srv_f(i);
    totn_srv1(2,i) = natlength(2,i)*srv_q_m(i)*sel_srv_m(i);
  }
  pred_bio.initialize();
  fspbio.initialize();
  mspbio.initialize(); 
  for (i=styr;i<=endyr;i++)
  {
    fspbio(i) = natlength_mat(1,i)*wtf(2);
    mspbio(i) = natlength_mat(2,i)*wtm;
    // Selection pattern
//    if (i<1974) sel_srv_use = sel_srv1;
    fspbio_srv1(i) = q1*natlength_mat(1,i)*elem_prod(wtf(2),sel_srv_f(i));
    mspbio_srv1(i) = q1*natlength_mat(2,i)*elem_prod(wtm,   sel_srv_m(i));

    // this is predicted survey in numbers not biomass-don't adjust by max selectivity 
    totn_srv1(1,i) = natlength(1,i)*sel_srv_f(i);
    totn_srv1(2,i) = natlength(2,i)*sel_srv_m(i);
    totSrvNum = totn_srv1(1,i) + totn_srv1(2,i)+1e-7;

    pred_p_srv1_len_new(1,1,i) = elem_prod(sel_srv_f(i),natlength_inew(1,i))/totSrvNum;
    pred_p_srv1_len_old(1,1,i) = elem_prod(sel_srv_f(i),natlength_iold(1,i))/totSrvNum;
    pred_p_srv1_len_new(2,1,i) = elem_prod(sel_srv_f(i),natlength_mnew(1,i))/totSrvNum;
    pred_p_srv1_len_old(2,1,i) = elem_prod(sel_srv_f(i),natlength_mold(1,i))/totSrvNum;
    pred_p_srv1_len_new(1,2,i) = elem_prod(sel_srv_m(i),natlength_inew(2,i))/totSrvNum;
    pred_p_srv1_len_old(1,2,i) = elem_prod(sel_srv_m(i),natlength_iold(2,i))/totSrvNum;
    pred_p_srv1_len_new(2,2,i) = elem_prod(sel_srv_m(i),natlength_mnew(2,i))/totSrvNum;
    pred_p_srv1_len_old(2,2,i) = elem_prod(sel_srv_m(i),natlength_mold(2,i))/totSrvNum;

    pred_bio(i) += natlength_inew(1,i)*wtf(1)+
                  (natlength_mnew(1,i)+natlength_mold(1,i))*wtf(2)+
                  (natlength_inew(2,i)+natlength_mnew(2,i)+natlength_mold(2,i))*wtm;
  }  
  depletion = pred_bio(endyr) / pred_bio(styr);
  fspbios = fspbio(1974,endyr);
  mspbios = mspbio(1974,endyr);
  // Legal males
  legal_males.initialize();
  legal_srv_males.initialize();
  for (i=styr;i<=endyr;i++)
  {
    // legal is >=138mm take half the numbers in the 135-139 bin
    legal_males(i)=0.5*natlength(2,i,23);
    legal_males(i) += sum(natlength(2,i)(24,nlenm));
    legal_srv_males(i)  = 0.5*natlength(2,i,23)*sel_srv_m(i,23);
    legal_srv_males(i) += sum(elem_prod(natlength(2,i)(24,nlenm),sel_srv_m(i)(24,nlenm)));
  }
  legal_malesd = legal_males(1974,endyr);
  rec_early_sd=mfexp(mean_log_rec1_early+ rec_devf_early);
  recf_sd = mfexp(mean_log_rec(1)+rec_dev(1)(1974,endyr-1));
  recm_sd = mfexp(mean_log_rec(2)+rec_dev(2)(1974,endyr-1));
// ---------------------------------------------------------------------------

FUNCTION get_catch_at_len
  int sex,i,j;
  dvariable totn_trawl,popn_disc,popn_lmale,popn_fit;
  dvariable ratio1, ratio2;

  pred_catch.initialize();
  pred_catch_ret.initialize();
  pred_catch_snowd.initialize();
  pred_catch_rkd.initialize();
  pred_catch_female_d.initialize();
  pred_catch_female_snowd.initialize();
  pred_catch_female_rkd.initialize();
  pred_catch_trawl.initialize();
//   cout<<" to get catch at length "<<endl;
  for (i=styr;i<=endyr;i++)
  {
    for (j=1;j<= nlenm;j++)
    {
    //total male directed catch
      ratio1 = (F(1,i,j)/(F(1,i,j)+Fdisct(2,i,j)+Fdisc_snow(2,i,j)+Fdisc_rk(2,i,j))) * (1-Smat(2,1,i,j));
      catch_lmale_new(i,j) = ratio1*(natl_inew_fishtime(2,i,j) + natl_mnew_fishtime(2,i,j));
      ratio2 = F(2,i,j)/(F(2,i,j)+Fdisct(2,i,j)+Fdisc_snow(2,i,j)+Fdisc_rk(2,i,j)) * (1-Smat(2,2,i,j));
      catch_lmale_old(i,j) = ratio2*(natl_iold_fishtime(2,i,j) + natl_mold_fishtime(2,i,j));
      catch_lmale(i,j) = catch_lmale_new(i,j)+catch_lmale_old(i,j);
      pred_catch(i) += catch_lmale(i,j)*wtm(j);
      //retained male directed catch     
      ratio1 = F_ret(1,i,j)/(F(1,i,j)+Fdisct(2,i,j)+Fdisc_snow(2,i,j)+Fdisc_rk(2,i,j)) * (1-Smat(2,1,i,j));
      catch_male_ret_new(i,j) = ratio1*(natl_inew_fishtime(2,i,j) + natl_mnew_fishtime(2,i,j)); 
      ratio2 = F_ret(2,i,j)/(F(2,i,j)+Fdisct(2,i,j)+Fdisc_snow(2,i,j)+Fdisc_rk(2,i,j)) * (1-Smat(2,2,i,j));
      catch_male_ret_old(i,j) = ratio2*(natl_iold_fishtime(2,i,j) + natl_mold_fishtime(2,i,j));
      catch_male_ret(i,j) = catch_male_ret_new(i,j)+catch_male_ret_old(i,j);
      pred_catch_ret(i) += catch_male_ret(i,j)*wtm(j);
      //snow crab discard catch male
      ratio1 = Fdisc_snow(2,i,j)/(F(1,i,j)+Fdisct(2,i,j)+Fdisc_snow(2,i,j)+Fdisc_rk(2,i,j)) * (1-Smat(2,1,i,j));
      catch_male_snowd_new(i,j) = ratio1*(natl_inew_fishtime(2,i,j) + natl_mnew_fishtime(2,i,j)); 
      ratio2 = Fdisc_snow(2,i,j)/(F(2,i,j)+Fdisct(2,i,j)+Fdisc_snow(2,i,j)+Fdisc_rk(2,i,j)) * (1-Smat(2,2,i,j));
      catch_male_snowd_old(i,j) = ratio2*(natl_iold_fishtime(2,i,j) + natl_mold_fishtime(2,i,j));
      catch_male_snowd(i,j)= catch_male_snowd_new(i,j)+catch_male_snowd_old(i,j);
      pred_catch_snowd(i) += catch_male_snowd(i,j)*wtm(j);
      //red king crab discard catch male     
      ratio1 = Fdisc_rk(2,i,j)/(F(1,i,j)+Fdisct(2,i,j)+Fdisc_snow(2,i,j)+Fdisc_rk(2,i,j)) * (1-Smat(2,1,i,j));
      catch_male_rkd_new(i,j) = ratio1*(natl_inew_fishtime(2,i,j) + natl_mnew_fishtime(2,i,j)); 
      ratio2 = Fdisc_rk(2,i,j)/(F(2,i,j)+Fdisct(2,i,j)+Fdisc_snow(2,i,j)+Fdisc_rk(2,i,j)) * (1-Smat(2,2,i,j));
      catch_male_rkd_old(i,j) = ratio2*(natl_iold_fishtime(2,i,j) + natl_mold_fishtime(2,i,j));
      catch_male_rkd(i,j) = catch_male_rkd_new(i,j)+catch_male_rkd_old(i,j);
      pred_catch_rkd(i) += catch_male_rkd(i,j)*wtm(j);
      //directed tanner discard catch female
      ratio1 = Fdiscf(i,j)/(Fdiscf(i,j)+Fdisct(1,i,j)+Fdisc_snow(1,i,j)+Fdisc_rk(1,i,j)) * (1-Smat(1,1,i,j));
      catch_female_d_new(i,j) = ratio1*(natl_inew_fishtime(1,i,j) + natl_mnew_fishtime(1,i,j)); 
      ratio2 = Fdiscf(i,j)/(Fdiscf(i,j)+Fdisct(1,i,j)+Fdisc_snow(1,i,j)+Fdisc_rk(1,i,j)) * (1-Smat(1,2,i,j));
      catch_female_d_old(i,j) = ratio2*(natl_iold_fishtime(1,i,j) + natl_mold_fishtime(1,i,j));
      catch_female_d(i,j) = catch_female_d_new(i,j)+catch_female_d_old(i,j);
      pred_catch_female_d(i) += catch_female_d(i,j)*wtf(2,j);
      //snow crab discard catch female
      ratio1 = Fdisc_snow(1,i,j)/(Fdiscf(i,j)+Fdisct(1,i,j)+Fdisc_snow(1,i,j)+Fdisc_rk(1,i,j)) * (1-Smat(1,1,i,j));
      catch_female_snowd_new(i,j) = ratio1*(natl_inew_fishtime(1,i,j) + natl_mnew_fishtime(1,i,j)); 
      ratio2 = Fdisc_snow(1,i,j)/(Fdiscf(i,j)+Fdisct(1,i,j)+Fdisc_snow(1,i,j)+Fdisc_rk(1,i,j)) * (1-Smat(1,2,i,j));
      catch_female_snowd_old(i,j) = ratio2*(natl_iold_fishtime(1,i,j) + natl_mold_fishtime(1,i,j));
      catch_female_snowd(i,j) = catch_female_snowd_new(i,j)+catch_female_snowd_old(i,j);
      pred_catch_female_snowd(i) += catch_female_snowd(i,j)*wtf(2,j);
      //red king crab discard catch female
      ratio1 = Fdisc_rk(1,i,j)/(Fdiscf(i,j)+Fdisct(1,i,j)+Fdisc_snow(1,i,j)+Fdisc_rk(1,i,j)) * (1-Smat(1,1,i,j));
      catch_female_rkd_new(i,j) = ratio1*(natl_inew_fishtime(1,i,j) + natl_mnew_fishtime(1,i,j)); 
      ratio2 = Fdisc_rk(1,i,j)/(Fdiscf(i,j)+Fdisct(1,i,j)+Fdisc_snow(1,i,j)+Fdisc_rk(1,i,j)) * (1-Smat(1,2,i,j));
      catch_female_rkd_old(i,j) = ratio2*(natl_iold_fishtime(1,i,j) + natl_mold_fishtime(1,i,j));
      catch_female_rkd(i,j) = catch_female_rkd_new(i,j)+catch_female_rkd_old(i,j);
      pred_catch_female_rkd(i) += catch_female_rkd(i,j)*wtf(2,j);
      //trawl bycatch
      ratio1 = Fdisct(1,i,j)/(Fdiscf(i,j)+Fdisct(1,i,j)+Fdisc_snow(1,i,j)+Fdisc_rk(1,i,j)) * (1-Smat(1,1,i,j));
      catch_trawl_female(i,j) = ratio1*(natl_inew_fishtime(1,i,j) + natl_mnew_fishtime(1,i,j)+natl_iold_fishtime(1,i,j) + natl_mold_fishtime(1,i,j)); 
      pred_catch_trawl(i) += catch_trawl_female(i,j)*wtf(2,j);
      ratio1 = Fdisct(2,i,j)/(F(1,i,j)+Fdisct(2,i,j)+Fdisc_snow(2,i,j)+Fdisc_rk(2,i,j)) * (1-Smat(2,1,i,j));
      catch_trawl_male(i,j) = ratio1*(natl_inew_fishtime(2,i,j) + natl_mnew_fishtime(2,i,j)+natl_iold_fishtime(2,i,j) + natl_mold_fishtime(2,i,j)); 
      pred_catch_trawl(i) += catch_trawl_male(i,j)*wtm(j);
    }   //end of J loop
  }    //end of year loop
  for (i=styr;i<=endyr;i++)
  { 
    // Discard catch by sex (AEP assumes all females are mature - for weight purposes)
    pred_catch_disc(1,i) = pred_catch_female_d(i);
    // Retained catch (males)
    popn_fit = sum(catch_male_ret(i)+1e-6);
    pred_p_fish_fit(1,i) = catch_male_ret_new(i)/popn_fit;
    pred_p_fish_fit(2,i) = catch_male_ret_old(i)/popn_fit;
    // Total catch (males)
    popn_lmale = sum(catch_lmale(i)+1e-6);
    pred_p_fish(1,i) = catch_lmale_new(i)/popn_lmale;
    pred_p_fish(2,i) = catch_lmale_old(i)/popn_lmale;
    // female discards
    pred_p_fish_discf(i) = catch_female_d(i)/sum(catch_female_d(i)+1e-6);
    // snow crab discards female male
    pred_p_snow(1,i) = catch_female_snowd(i)/sum(catch_female_snowd(i)+1e-6);
    pred_p_snow(2,i) = catch_male_snowd(i)  /sum(catch_male_snowd(i)+1e-6);
    // red king crab discards female male
    pred_p_rk(2,i) = catch_male_rkd(i)/sum(catch_male_rkd(i)+1e-6);
    pred_p_rk(1,i) = catch_female_rkd(i)/sum(catch_female_rkd(i)+1e-6);
    // total trawl selected numbers
    totn_trawl = sum(catch_trawl_male(i))+sum(catch_trawl_female(i)+1e-6);
    // Trawl proportions (adds to 1 over sex, shell and length)
    pred_p_trawl(1,i) = catch_trawl_female(i)/totn_trawl;
    pred_p_trawl(2,i) = catch_trawl_male(i)  /totn_trawl;
  }
// ---------------------------------------------------------------------------
                                                               
FUNCTION evaluate_the_objective_function
  int yrc, ii, ij, mat, sex, ik, i, j, k;
  dvariable nextf;
  dvar_vector tmpsum(1,nlenm);
  dvar_vector tmpsum2(1,nlenm);
  dvariable multi;
  like_initn=0.0;
  tmpsum.initialize();
  tmpsum2.initialize();
  Fout.initialize();
  // PENALTIES
  // =========
  // Constrains on recruitment
  penal_rec.initialize();
  if (active(rec_devf))
  {
    //recruitment likelihood - norm2 is sum of square values   
    penal_rec  = 1.0*like_wght_recf*norm2(rec_devf); //+ like_wght_rec*norm2(rec_devm);
    penal_rec += 1.0*norm2(first_difference(rec_devf_early));
    Fout(1) = penal_rec; f += penal_rec; 
    // Deviations ondifference between make and female recruites (not used)
    penal_sexr.initialize();
    if (active(rec_devm))
    {
      for(i=styr;i<endyr;i++)
        penal_sexr += like_wght_sexr*square((mean_log_rec(1)+rec_devf(i))-(mean_log_rec(2)+rec_devm(i)));
    }
    Fout(2) = penal_sexr; f += penal_sexr;
  } 
  //nat Mort. penalty
  if(active(Mmult_imat))
  {  
    nat_penalty = 0.5 * square((Mmult_imat - 1.0) / 0.05);
    f += nat_penalty; Fout(27) += nat_penalty;
  }
  if(active(Mmultm))
  {  
    nat_penalty = 0.5 * square((Mmultm - 1.0) / 0.05);
    f += nat_penalty; Fout(27) += nat_penalty;
  }
  if(active(Mmultf))
  {  
    nat_penalty = 0.5 * square((Mmultf - 1.0) / 0.05);
    f += nat_penalty; Fout(27) += nat_penalty;
  }

//penalty on survey Q
  if (q_prior_switch==1)
  {
    // if(active(srv3_q))
    // {  
      //max of underbag at 182.5 mm is 0.873   
      // srv3q_penalty = 0.5 * square((srv3_q - 0.88) / 0.05);
      // f += srv3q_penalty; 
      // Fout(4) = srv3q_penalty;
    // }
    // if(active(srv3_femQ))
    // {  
    //peak of females is at about 80mm underbag is 0.75 at this size - less uncertainty  
      // srv3q_penalty = 0.5 * square((srv3_femQ - 0.88) / 0.05);
      // f += srv3q_penalty; 
      // Fout(4) += srv3q_penalty;
    // }
  }
  
  // Bayesian part - likelihood on growth parameters af,am,bf,bm
  // not used in this case
  af_penal = 0; bf_penal = 0; am_penal = 0; bm_penal = 0;
  if(active(af1))
  {  
    bf_penal = 0.5 * square((bf1 - 0.9132661) / 0.025);
    f += bf_penal; Fout(5) = bf_penal;
    af_penal = 0.5 * square((af1 - 0.56560241)    / 0.1);
    f += af_penal; Fout(5) += af_penal;
  }
  if(active(am1))
  {
    am_penal   = 0.5 * square((am1 - 0.437941)/0.025);
    f += am_penal; Fout(6) = am_penal;
  }
  if(active(bm1))
  {
    //am_penal   = 0.5 * square((am1 - 0.437941)/0.025);
    bm_penal = 0.5 * square((bm1 - 0.9487) /0.1);
    f += bm_penal; Fout(7) = bm_penal;
  }
  if(active(matestm))
  {
    like_mat = 1.0*norm2(first_difference(first_difference(matestf)));
    like_mat += 0.5*norm2(first_difference(first_difference(matestm)));
//    like_mat += 2.0*norm2(maturity_est(2)-maturity_logistic)/(0.05*0.05);
//    like_mat += 1.0*norm2(maturity_est(1)-maturity_average(1))/(0.05*0.05);
    f += like_mat;
    Fout(3)= like_mat;
  }
  
  sel_50m_penal = 0;
  if(active(log_sel50_dev_mn))
    sel_50m_penal = like_wght_sel50*norm2(first_difference(log_sel50_dev_mn));
  Fout(9) = sel_50m_penal; f += sel_50m_penal; 

  // penalty on smoothness of initial year length comp
//  initsmo_penal = norm2(first_difference(mnatlen_styr(1)))+norm2(first_difference(mnatlen_styr(2)))+
//                 norm2(first_difference(fnatlen_styr(1)))+norm2(first_difference(fnatlen_styr(2)));
//  initsmo_penal = norm2(first_difference(mnatlen_styr))+norm2(first_difference(fnatlen_styr));
//  Fout(10) = initsmo_penal; f += initsmo_penal;
  
  // various penalties
  // =================
  fmort_pen.initialize();
  if (!last_phase())
  {
    if(active(fmort_dev))       fmort_pen(1) += 10. *norm2(fmort_dev);
    if(active(fmortd_snow_dev)) fmort_pen(2) += 5.0 *norm2(fmortd_snow_dev);
    if(active(fmortt_dev))      fmort_pen(3) += 5.0 *norm2(fmortt_dev);
    if(active(fmortd_rk_dev))   fmort_pen(4) += 3.0 *norm2(fmortd_rk_dev);
  }
  else
  {
    if(active(fmort_dev))       fmort_pen(1) += 1.0 *norm2(fmort_dev);
    if(active(fmortd_snow_dev)) fmort_pen(2) += 1.0 *norm2(fmortd_snow_dev);
    if(active(fmortt_dev))      fmort_pen(3) += 1.0 *norm2(fmortt_dev);
    if(active(fmortd_rk_dev))   fmort_pen(4) += 1.0 *norm2(fmortd_rk_dev);
  }

  fpen = sum(fmort_pen);
  f += fpen; Fout(11) = fpen;
  
  // LIKELIHOODS
  // ===========
  
  len_like.initialize();
  f += fpen; Fout(11) = fpen;
  
  // LIKELIHOODS
  // ===========
  len_like.initialize();
  // fishery length likelihood (old and new shell together)
  for (i=1; i <= nobs_fish; i++)
  {
    ii = yrs_fish(i);
    for (j=1; j<=nlenm; j++){
      len_like(1) -= nsamples_fish(1,i)*(obs_p_fish_ret(1,i,j)+obs_p_fish_ret(2,i,j))*
                    log(pred_p_fish_fit(2,ii,j)+pred_p_fish_fit(1,ii,j)+p_const);
    }
    tmpsum+=obs_p_fish_ret(1,i)+obs_p_fish_ret(2,i);
    tmpsum2+=pred_p_fish_fit(2,ii)+pred_p_fish_fit(1,ii);
  }
  for (i=1; i <= nobs_fish_discm; i++)
  {
    ij = yrs_fish_discm(i);
    for (j=1; j<=nlenm; j++){
      len_like(2) -= nsamples_fish_discm(1,i)*(obs_p_fish_tot(1,i,j)+obs_p_fish_tot(2,i,j))*log(pred_p_fish(1,ij,j)+pred_p_fish(2,ij,j)+p_const);
    }
  }
  // fishery length likelihood (female discards)
  for (i=1; i <= nobs_fish_discf; i++)
  {
    ik=yrs_fish_discf(i);
    for (j=1; j<=nlenm; j++){
      len_like(3) -= nsamples_fish_discf(i)*obs_p_fish_discf(i,j)*log(pred_p_fish_discf(ik,j)+p_const);
    }
  }
  // trawl fishery length likelihood 
  for (i=1; i <= nobs_trawl; i++)
  {
    ij=yrs_trawl(i);
    for(sex=1;sex<=2;sex++){
      for (j=1; j<=nlenm; j++){
        len_like(5) -= nsamples_trawl(sex,i)*obs_p_trawl(sex,i,j)*log(pred_p_trawl(sex,ij,j)+p_const);
   }}}
  // fishery length likelihood snow crab fishery discards
  for (i=1; i <= nobs_snowfish_discf; i++)
  {
    ik=yrs_snowfish_discf(i);
    for (j=1; j<=nlenm; j++){
      len_like(6) -= nsamples_snowfish_discf(i)*obs_p_snow(1,i,j)*log(pred_p_snow(1,ik,j)+p_const);
    }
  }
   
  for (i=1; i <= nobs_snowfish_discm; i++)
  {
    ik=yrs_snowfish_discm(i);
    for (j=1; j<=nlenm; j++){
      len_like(6) -= nsamples_snowfish_discm(1,i)*obs_p_snow(2,i,j)*log(pred_p_snow(2,ik,j)+p_const);
    }
  }
  // fishery length likelihood red king crab discards
  for (i=1; i <= nobs_rkfish_discf; i++)
  {
    ik = yrs_rkfish_discf(i);
    for (j=1; j<=nlenm; j++)
    {
      len_like(7) -= nsamples_rkfish_discf(i)  *obs_p_rk(1,i,j)*log(pred_p_rk(1,ik,j)+p_const);
      len_like(8) -= nsamples_rkfish_discm(1,i)*obs_p_rk(2,i,j)*log(pred_p_rk(2,ik,j)+p_const);
    }
  }
  
  //add the offset to the likelihood   
  len_like(1) -= offset(1);
  len_like(2) -= offset(2);
  len_like(3) -= offset(3);
  len_like(5) -= offset(5);
  len_like(6) -= offset(6);
  len_like(7) -= offset(7);
  len_like(8) -= offset(8);
  
  // survey likelihood
  for(k=1;k<=2;k++){  //sex
    for (i=1; i <=nobs_srv1_length; i++)
    {
      ii=yrs_srv1_length(i);
      for (j=1; j<=nlenm; j++)
      {
        // obs(maturity, SC, sex, year), pred(maturity,sex, year)
        // this is for mature new and old shell together
        len_like(4) -= nsamples_srv1_length(2,1,k,i)*(
                       obs_p_srv1_len(2,1,k,i,j)+
                       obs_p_srv1_len(2,2,k,i,j))*
                       log(pred_p_srv1_len_new(2,k,ii,j)+
                           pred_p_srv1_len_old(2,k,ii,j)+p_const);
        // immature new and old together
        len_like(4) -= nsamples_srv1_length(2,2,k,i)*(obs_p_srv1_len(1,1,k,i,j)+obs_p_srv1_len(1,2,k,i,j))*log(pred_p_srv1_len_new(1,k,ii,j)+pred_p_srv1_len_old(1,k,ii,j)+p_const);
       }  //j loop     
     }// year loop
   } //sex loop
   len_like(4) -= offset(4);
  
//  if (current_phase() > 6) //?? AEP (weights are getting messy here - this should ALL be done in the input files)
  Fout(12) = like_wght(1)*len_like(1); // Retained fishery 
  Fout(13) = like_wght(2)*len_like(2); // Total (ret+disc)
  Fout(14) = like_wght(3)*len_like(3); // Female
  Fout(15) = like_wght(4)*len_like(4); // Survey
  Fout(16) = like_wght(7)*len_like(5); // Trawl
  Fout(17) =              len_like(6); // snow fishery
  Fout(18) = len_like(7) +len_like(8);// red king fishery male plus female
  f += sum(Fout(12,18));

  // this fits mature biomass separate male and female
  //female biomass only for 1974 to endyr, male biomass from 1969 to endyr
  surv_like  = norm2(elem_div( log(obs_srv1_spbiom(1)(yrs_srv1)+.000001)-log(fspbio_srv1(yrs_srv1)+.000001),
               sqrt(2)*sqrt(log(elem_prod(cv_srv1(1)(yrs_srv1),cv_srv1(1)(yrs_srv1))+1.0))));
  surv_like += norm2(elem_div( log(obs_srv1_spbiom(2)(yrs_srv1)+.000001)-log(mspbio_srv1(yrs_srv1)+.000001),
               sqrt(2)*sqrt(log(elem_prod(cv_srv1(2)(yrs_srv1),cv_srv1(2)(yrs_srv1))+1.0))));

  Fout(19) = surv_like;  
  f       += surv_like; 

  //fishery cpue likelihood
  cpue_pred = cpueq * legal_males;
  cpue_like=0.0;
  for(yrc=1969;yrc<=endyr-1;yrc++)
  {
    // don't use cpue in years when fishery was closed
    if(cpue(yrc+1)>0.0001)
    {
      cpue_like += 0.001*square( log((cpue(yrc+1)+1e-6)/(cpue_pred(yrc)+1e-6)) / cpue_cv);
    }
  }
  // ===============================================================
  // catch likelihoods
  // don't include last year as that would be endyr+1 fishery season
  // ===============================================================
  catch_like1=.0;
  catch_like2=.0;
  catch_likef=.0;
  catch_liket=.0;
  catch_likes=.0;
  catch_liker=.0;

  catch_like1 = norm2(log(obs_catchtot_biom(yrs_fish_catchf)+1e-5)-log(pred_catch(yrs_fish_catchf)        +1e-5));
  catch_likef = norm2(log(obs_catchdf_biom(yrs_fish_catchf) +1e-5)-log(pred_catch_disc(1)(yrs_fish_catchf)+1e-5));
  catch_like2 = norm2(log(catch_ret(1965,endyr-1)+1e-5)           -log(pred_catch_ret(1965,endyr-1)       +1e-5));
  catch_liket = norm2(log(obs_catcht_biom(yrs_trawl_c)+1e-5)      -log(pred_catch_trawl(yrs_trawl_c)      +1e-5));
  for(i=1;i<=nobs_discardc;i++){
    pred_tmp(1,i)=pred_catch_female_snowd(yrs_discardc(i));
    pred_tmp(2,i)=pred_catch_snowd(yrs_discardc(i));
    pred_tmp(3,i)=pred_catch_rkd(yrs_discardc(i));
    pred_tmp(4,i)=pred_catch_female_rkd(yrs_discardc(i));
  }
  catch_likes  = norm2(log(catch_snowodisc(1)+1e-5)-log(pred_tmp(1)+1e-5));
  catch_likes += norm2(log(catch_snowodisc(2)+1e-5)-log(pred_tmp(2)+1e-5));
  catch_liker  = norm2(log(catch_rkodisc(2)  +1e-5)-log(pred_tmp(3)+1e-5));
  catch_liker += norm2(log(catch_rkodisc(1)  +1e-5)-log(pred_tmp(4)+1e-5));
  Fout(21) = like_wght(6)*catch_like1;
  Fout(22) = like_wght(6)*catch_like2;
  Fout(23) = like_wght(6)*catch_likef;
  Fout(24) = like_wght(6)*catch_liket;
  Fout(25) = like_wght(6)*catch_likes;
  Fout(26) = like_wght(6)*catch_liker;
  f += sum(Fout(21,26));
  call_no += 1;

// ========================y==================================================

FUNCTION Misc_output
  int i,j,ii,ij,k,sex;
  dvariable tmpi;
  dvar_matrix sel_srv_use(1,2,1,nlenm);
  dvar_matrix cv_srv1_nowt(1,2,styr,endyr);
//legal size for tanner is 138mm
  pred_catch_gt101.initialize();
  pred_catch_no_gt101.initialize();

//   cout<<" to misc output "<<endl;
  //weight each years estimate by 1/(2*variance) - use cv of biomass in sqrt(log(cv^2+1)) as sd of log(biomass) 
  for(i=1;i<=nobs_srv1;i++)
  {
    cv_srv1_nowt(1,yrs_srv1(i)) = cv_srv1o(1,i);
    cv_srv1_nowt(2,yrs_srv1(i)) = cv_srv1o(2,i);
    biom_tmp(1,yrs_srv1(i)) = fspbio_srv1(yrs_srv1(i));
    biom_tmp(2,yrs_srv1(i)) = mspbio_srv1(yrs_srv1(i));
  }
  
  // Combined likelihood
  surv_like_nowt.initialize();
  for(sex=1;sex<=2;sex++)
  {
    surv_like_nowt += norm2(elem_div( log(obs_srv1_bioms(sex)(yrs_srv1)+.01)-log(biom_tmp(sex)(yrs_srv1)+.01),
                sqrt(2)*sqrt(log(elem_prod(cv_srv1_nowt(sex)(yrs_srv1),cv_srv1_nowt(sex)(yrs_srv1))+1.0))));
  }
  for (i=styr;i<=endyr;i++)
    for (j = 23 ; j<= nlenm; j++)
    {
       pred_catch_no_gt101(i) += (F(1,i,j)/(F(1,i,j)+Fdisct(2,i,j)))*natl_inew_fishtime(2,i,j)*(1-Smat(2,1,i,j)) + 
                               (F(1,i,j)/(F(1,i,j)+Fdisct(2,i,j)))*natl_mnew_fishtime(2,i,j)*(1-Smat(2,1,i,j))+ 
                               (F(2,i,j)/(F(2,i,j)+Fdisct(2,i,j)))*natl_iold_fishtime(2,i,j)*(1-Smat(2,2,i,j))+
                               (F(2,i,j)/(F(2,i,j)+Fdisct(2,i,j)))*natl_mold_fishtime(2,i,j)*(1-Smat(2,2,i,j));
       pred_catch_gt101(i)+= (F(1,i,j)/(F(1,i,j)+Fdisct(2,i,j)))*natl_inew_fishtime(2,i,j)*(1-Smat(2,1,i,j)) + 
                           (F(1,i,j)/(F(1,i,j)+Fdisct(2,i,j)))*natl_mnew_fishtime(2,i,j)*(1-Smat(2,1,i,j))+ 
                           (F(2,i,j)/(F(2,i,j)+Fdisct(2,i,j)))*natl_iold_fishtime(2,i,j)*(1-Smat(2,2,i,j))+
                           (F(2,i,j)/(F(2,i,j)+Fdisct(2,i,j)))*natl_mold_fishtime(2,i,j)*(1-Smat(2,2,i,j)) * wtm(j);
     if (j<24) // AEP???
       {
        pred_catch_gt101(i)=pred_catch_gt101(i)*0.5;
        pred_catch_no_gt101(i)=pred_catch_no_gt101(i)*0.5;
       }
    }
// cout<<" to large males "<<endl;
  bio_males_gt101.initialize();
  num_males_gt101.initialize();
  for (i=styr;i<=endyr;i++)
   for(j=23;j<=nlenm;j++)
    {
     num_males_gt101(i)+= natl_inew_fishtime(2,i,j) + natl_iold_fishtime(2,i,j) + natl_mnew_fishtime(2,i,j) + natl_mold_fishtime(2,i,j);
     bio_males_gt101(i)+= (natl_inew_fishtime(2,i,j) + natl_iold_fishtime(2,i,j) + natl_mnew_fishtime(2,i,j) + natl_mold_fishtime(2,i,j))*wtm(j);
     if (j<24)
      {
       num_males_gt101(i)=num_males_gt101(i)*0.5;
       bio_males_gt101(i)=bio_males_gt101(i)*0.5;
      }
    }
//  cout<<" to eff N "<<endl;
  // Effective N's
  for (i=1; i <= nobs_fish; i++)
   for(k=1;k<=2;k++)
    {
     ii=yrs_fish(i);
    if(sum(obs_p_fish_ret(k,i))<0.00001)
     effn_fish_ret(k,ii)=1/(norm2(pred_p_fish_fit(k,ii)-obs_p_fish_ret(k,i))/(pred_p_fish_fit(1,ii)*(1-pred_p_fish_fit(1,ii))+pred_p_fish_fit(2,ii)*(1-pred_p_fish_fit(2,ii))));
    } 
  for (i=1; i <= nobs_fish_discm; i++)
   {
    ij=yrs_fish_discm(i);
    for(k=1;k<=2;k++)
        if(sum(obs_p_fish_tot(k,i))<0.00001)
     effn_fish_tot(k,ij)=1/(norm2(pred_p_fish(k,ij)-obs_p_fish_tot(k,i))/(pred_p_fish(1,ij)*(1-pred_p_fish(1,ij))+pred_p_fish(2,ij)*(1-pred_p_fish(2,ij))));
   }
  for(k=1;k<=2;k++)  //sex
   for (i=1; i <=nobs_srv1_length; i++)
    {
     ii=yrs_srv1_length(i);
     effn_srv1(1,1,k,ii)=1./(norm2(pred_p_srv1_len_new(1,k,ii)-obs_p_srv1_len(1,1,k,i))/(pred_p_srv1_len_new(1,1,ii)*(1-pred_p_srv1_len_new(1,1,ii))+pred_p_srv1_len_new(1,2,ii)*(1-pred_p_srv1_len_new(1,2,ii))+pred_p_srv1_len_old(1,1,ii)*(1-pred_p_srv1_len_old(1,1,ii))+pred_p_srv1_len_old(1,2,ii)*(1-pred_p_srv1_len_old(1,2,ii))));
     if(k > 1) 
       effn_srv1(1,2,k,ii)=0.0;
     else
       effn_srv1(1,2,k,ii)=0.0;
     effn_srv1(2,1,k,ii)=1./(norm2(pred_p_srv1_len_new(2,k,ii)-(obs_p_srv1_len(2,1,k,i)))/(pred_p_srv1_len_new(2,1,ii)*(1-pred_p_srv1_len_new(2,1,ii))+pred_p_srv1_len_new(2,2,ii)*(1-pred_p_srv1_len_new(2,2,ii))+pred_p_srv1_len_old(2,1,ii)*(1-pred_p_srv1_len_old(2,1,ii))+pred_p_srv1_len_old(2,2,ii)*(1-pred_p_srv1_len_old(2,2,ii))));
     effn_srv1(2,2,k,ii)=1./(norm2(pred_p_srv1_len_old(2,k,ii)-(obs_p_srv1_len(2,2,k,i)))/(pred_p_srv1_len_new(2,1,ii)*(1-pred_p_srv1_len_new(2,1,ii))+pred_p_srv1_len_new(2,2,ii)*(1-pred_p_srv1_len_new(2,2,ii))+pred_p_srv1_len_old(2,1,ii)*(1-pred_p_srv1_len_old(2,1,ii))+pred_p_srv1_len_old(2,2,ii)*(1-pred_p_srv1_len_old(2,2,ii))));
 
    } // year loop

  // spawning biomass and related outputs
  efspbio_matetime.initialize();
  emspbio_matetime.initialize();
  mspbio_old_matetime.initialize();
  fspbio_new_matetime.initialize();
  efspbio_new_matetime.initialize();
  fspnum_new_matetime.initialize();
  efspnum_matetime.initialize();
  emspnum_old_matetime.initialize();
  mspnum_matetime.initialize();
  for (i=styr;i<= endyr;i++)
  {
    mspbio_matetime(i) = (elem_prod(Smat(2,1,i)*mfexp(-spmo*M_matn(2)),mfexp(-catch_midpt(i)*M_matn(2))*natlength_mnew(2,i))+elem_prod(Smat(2,2,i)*mfexp(-spmo*M_mato(2)),mfexp(-catch_midpt(i)*M_mato(2))*natlength_mold(2,i)))*wtm;
    fspbio_matetime(i) = (elem_prod(Smat(1,1,i)*mfexp(-spmo*M_matn(1)),mfexp(-catch_midpt(i)*M_matn(1))*natlength_mnew(1,i))+elem_prod(Smat(1,2,i)*mfexp(-spmo*M_mato(1)),mfexp(-catch_midpt(i)*M_mato(1))*natlength_mold(1,i)))*wtf(2);
    // Extra mortality
    if(i==1983 && killem_switch==1)
    {
      mspbio_matetime(i) = (elem_prod(Smat(2,1,i)*mfexp(-spmo*M_matn(2)*mat_big(2)),mfexp(-catch_midpt(i)*M_matn(2)*mat_big(2))*natlength_mnew(2,i))+elem_prod(Smat(2,2,i)*mfexp(-spmo*M_mato(2)*mat_big(2)),mfexp(-catch_midpt(i)*M_mato(2)*mat_big(2))*natlength_mold(2,i)))*wtm;
      fspbio_matetime(i) = (elem_prod(Smat(1,1,i)*mfexp(-spmo*M_matn(1)*mat_big(1)),mfexp(-catch_midpt(i)*M_matn(1)*mat_big(1))*natlength_mnew(1,i))+elem_prod(Smat(1,2,i)*mfexp(-spmo*M_mato(1)*mat_big(1)),mfexp(-catch_midpt(i)*M_mato(1)*mat_big(1))*natlength_mold(1,i)))*wtf(2);
    }
    else
    {
      mspbio_fishtime(i) = (natl_mnew_fishtime(2,i)+natl_mold_fishtime(2,i))*wtm;
      fspbio_fishtime(i) = (natl_mnew_fishtime(1,i)+natl_mold_fishtime(1,i))*wtf(2);
    }
    
    if(i==1983 && killem_switch==1)
    {
      for(j=1;j<=nlenm;j++)
      {
        emspnum_old_matetime(i) += Smat(2,2,i,j)*mfexp(-spmo*M_mato(2)*mat_big(2))*mfexp(-catch_midpt(i)*M_mato(2)*mat_big(2))*natlength_mold(2,i,j);
        mspnum_matetime(i) += Smat(2,1,i,j)*mfexp(-spmo*M_matn(2)*mat_big(2))*mfexp(-catch_midpt(i)*M_matn(2)*mat_big(2))*natlength_mnew(2,i,j) + 
                              Smat(2,2,i,j)*mfexp(-spmo*M_mato(2)*mat_big(2))*mfexp(-catch_midpt(i)*M_mato(2)*mat_big(2))*natlength_mold(2,i,j);
        mspbio_old_matetime(i) += (Smat(2,2,i,j)*mfexp(-spmo*M_mato(2)*mat_big(2))*mfexp(-catch_midpt(i)*M_mato(2)*mat_big(2))*natlength_mold(2,i,j))*wtm(j);
        fspnum_new_matetime(i) += (Smat(1,1,i,j)*mfexp(-spmo*M_matn(1)*mat_big(1))*mfexp(-catch_midpt(i)*M_matn(1)*mat_big(1))*natlength_mnew(1,i,j));
        fspbio_new_matetime(i)+= (Smat(1,1,i,j)*mfexp(-spmo*M_matn(1)*mat_big(1))*mfexp(-catch_midpt(i)*M_matn(1)*mat_big(1))*natlength_mnew(1,i,j))*wtf(2,j);
        efspnum_matetime(i) += (Smat(1,1,i,j)*mfexp(-spmo*M_matn(1)*mat_big(1))*mfexp(-catch_midpt(i)*M_matn(1)*mat_big(1))*natlength_mnew(1,i,j)+Smat(1,2,i,j)*mfexp(-spmo*M_mato(1)*mat_big(1))*mfexp(-catch_midpt(i)*M_mato(1)*mat_big(1))*natlength_mold(1,i,j));
      }
    }
    else
    {
      for(j=1;j<=nlenm;j++)
      {
        emspnum_old_matetime(i) += Smat(2,2,i,j)*mfexp(-spmo*M_mato(2))*mfexp(-catch_midpt(i)*M_mato(2))*natlength_mold(2,i,j);
        mspnum_matetime(i)      += Smat(2,1,i,j)*mfexp(-spmo*M_matn(2))*mfexp(-catch_midpt(i)*M_matn(2))*natlength_mnew(2,i,j) + 
                                   Smat(2,2,i,j)*mfexp(-spmo*M_mato(2))*mfexp(-catch_midpt(i)*M_mato(2))*natlength_mold(2,i,j);
        mspbio_old_matetime(i) += (Smat(2,2,i,j)*mfexp(-spmo*M_mato(2))*mfexp(-catch_midpt(i)*M_mato(2))*natlength_mold(2,i,j))*wtm(j);
        fspnum_new_matetime(i) += (Smat(1,1,i,j)*mfexp(-spmo*M_matn(1))*mfexp(-catch_midpt(i)*M_matn(1))*natlength_mnew(1,i,j));
        fspbio_new_matetime(i) += (Smat(1,1,i,j)*mfexp(-spmo*M_matn(1))*mfexp(-catch_midpt(i)*M_matn(1))*natlength_mnew(1,i,j))*wtf(2,j);
        efspnum_matetime(i)    += (Smat(1,1,i,j)*mfexp(-spmo*M_matn(1))*mfexp(-catch_midpt(i)*M_matn(1))*natlength_mnew(1,i,j)+Smat(1,2,i,j)*mfexp(-spmo*M_mato(1))*mfexp(-catch_midpt(i)*M_mato(1))*natlength_mold(1,i,j));
      }
    }
    // effective sp numbers
    emspbio_matetime(i) = mspbio_old_matetime(i);
      
    // for male old shell mating only (AEP ERROR IN OLD CODE HAS >=)
    efspbio_matetime(i) = fspbio_matetime(i);
    if (emspnum_old_matetime(i) < (efspnum_matetime(i)/mate_ratio))
      efspbio_matetime(i) = fspbio_matetime(i)*((emspnum_old_matetime(i)*mate_ratio)/efspnum_matetime(i));
       
    // effective sp numbers for new shell females
    efspbio_new_matetime(i) = fspbio_new_matetime(i);
    if (emspnum_old_matetime(i) < fspnum_new_matetime(i)/mate_ratio)
      efspbio_new_matetime(i) = fspbio_new_matetime(i)*((emspnum_old_matetime(i)*mate_ratio)/fspnum_new_matetime(i));
   }
  // Sex ratio
  for (i=styr;i<=endyr;i++)
    if((sum(natlength(1,i))+sum(natlength(2,i)))<0.01)
    { 
      predpop_sexr(i)=0.0;
    }
    else
    {
      predpop_sexr(i)=sum(natlength(1,i))/(sum(natlength(1,i))+sum(natlength(2,i)));
    }
  // Age-structure
  natlength_mold_age.initialize();
  
  // initialize
  tmpi = 1.0;
  for(j=1;j<=(nages-3);j++) 
    tmpi += mfexp(-j*M_mato(1));
  natlength_mold_age(1,styr,1) = natlength_mold(1,styr)/(tmpi+(mfexp(-(nages-2)*M_mato(1))/(1-mfexp(-M_mato(1)))));

  for(j=1;j<=(nages-2);j++)
    natlength_mold_age(1,styr,j+1) = natlength_mold_age(1,styr,1)*mfexp(-j*M_mato(1));

  natlength_mold_age(1,styr,nages) = natlength_mold_age(1,styr,1)*(mfexp(-(nages-2)*M_mato(1))/(1-mfexp(-M_mato(1))));

  tmpi = 1.0;
  for(j=1;j<=(nages-3);j++) 
    tmpi += mfexp(-j*M_mato(2));
  natlength_mold_age(2,styr,1) = natlength_mold(2,styr)/(tmpi+(mfexp(-(nages-2)*M_mato(2))/(1-mfexp(-M_mato(2)))));
  for(j=1;j<=(nages-2);j++)
    natlength_mold_age(2,styr,j+1) = natlength_mold_age(2,styr,1)*mfexp(-j*M_mato(2));
  natlength_mold_age(2,styr,nages) = natlength_mold_age(2,styr,1)*(mfexp(-(nages-2)*M_mato(2))/(1-mfexp(-M_mato(2))));

  //numbers at length from styr to endyr
  for(sex=1;sex<=2;sex++)
  {
    for (i=styr;i< endyr;i++)
    {
     // for numbers by length and age assumes no molting after maturity
      natlength_mold_age(sex,i+1,1) = mfexp(-(1-catch_midpt(i))*M_matn(sex)) * 
                                       elem_prod(Smat(sex,1,i),mfexp(-catch_midpt(i)*M_matn(sex))*natlength_mnew(sex,i));
      for(j=1;j<=(nages-1);j++)
        natlength_mold_age(sex,i+1,j+1) = (mfexp(-(1-catch_midpt(i))*M_mato(sex)) * 
                                           elem_prod(Smat(sex,2,i),mfexp(-catch_midpt(i)*M_mato(sex))*
                                           natlength_mold_age(sex,i,j)));

      natlength_mold_age(sex,i+1,nages) += (mfexp(-(1-catch_midpt(i))*M_mato(sex)) * 
                                   elem_prod(Smat(sex,2,i),mfexp(-catch_midpt(i)*M_mato(sex))*natlength_mold_age(sex,i,nages)));
                                   
    }
  }
  // Legal males
  popn.initialize();
  legal_srv_males_n.initialize();
  legal_srv_males_o.initialize();
  pred_srv1.initialize();
  pred_srv1_bioms.initialize();
  for (i=styr;i<=endyr;i++)
  {
    // Selection pattern
    // if (i<1978) sel_srv_use = sel_srv1;
    // legal is >102mm take half the numbers in the 100-105 bin
    legal_males_bio(i) = legal_males(i)*wtm(23);
    legal_srv_males_n(i) = 0.5*natlength_new(2,i,23)*sel_srv_m(i,23);
    legal_srv_males_o(i) = 0.5*natlength_old(2,i,23)*sel_srv_m(i,23);
    legal_srv_males_bio(i) = legal_srv_males(i)*wtm(23);
    for(j=24;j<=nlenm;j++)
    {
      legal_males_bio(i) += natlength(2,i,j)*wtm(j);
      legal_srv_males_n(i) += natlength_new(2,i,j)*sel_srv_use(1,j);
      legal_srv_males_o(i) += natlength_old(2,i,j)*sel_srv_use(1,j);
      legal_srv_males_bio(i) += natlength(2,i,j)*sel_srv_use(1,j)*wtm(j);
    }
    
    // survey numbers
    fspbio_srv1_num(1,i) = q1*natlength_mnew(1,i)*sel_srv_use(1);
    mspbio_srv1_num(1,i) = q1*natlength_mnew(2,i)*sel_srv_use(2);
    fspbio_srv1_num(2,i) = q1*natlength_mold(1,i)*sel_srv_use(1);
    mspbio_srv1_num(2,i) = q1*natlength_mold(2,i)*sel_srv_use(2);
    
    // total survey summaries
    for(sex=1;sex<=2;sex++)
    {
      if(sex<2)
        pred_srv1_bioms(sex,i) = q1*((natlength_inew(sex,i)*elem_prod(sel_srv_use(sex),wtf(1)))+
                                 ((natlength_mnew(sex,i)+natlength_mold(sex,i))*elem_prod(sel_srv_use(sex),wtf(2))));
      else
        pred_srv1_bioms(sex,i) = q1*(natlength(sex,i)*elem_prod(sel_srv_use(sex),wtm));

      pred_srv1(sex,i) = q1*elem_prod(natlength(sex,i),sel_srv_use(sex));
      popn(i) += sum(natlength(sex,i));
    } 
  }
  // Survey likelihood (by year)
  len_like_srv.initialize();
  for(sex=1;sex<=2;sex++) 
  {
    for (i=1; i <=nobs_srv1_length; i++)
    {
      ii=yrs_srv1_length(i);
     
      for (j=1; j<=nlenm; j++)
      {
       // immature new and old together in likelihood indices are (mat,shell,sex,year,length)
        len_like_srv(1,1,sex) -= nsamples_srv1_length(1,1,sex,i)*(1e-9+obs_p_srv1_len(1,1,sex,i,j)+obs_p_srv1_len(1,2,sex,i,j))*log(pred_p_srv1_len_new(1,sex,ii,j)+pred_p_srv1_len_old(1,sex,ii,j)+1e-9);
        len_like_srv(1,2,sex) = 0.0;
        
       // mature
        len_like_srv(2,1,sex) -= nsamples_srv1_length(2,1,sex,i)*(1e-9+obs_p_srv1_len(2,1,sex,i,j))*log(pred_p_srv1_len_new(2,sex,ii,j)+1e-9);
        len_like_srv(2,2,sex) -= nsamples_srv1_length(2,2,sex,i)*(1e-9+obs_p_srv1_len(2,2,sex,i,j))*log(pred_p_srv1_len_old(2,sex,ii,j)+1e-9);

      }  //j loop     
    } // year loop
  } // sex loop
// ==========================================================================

REPORT_SECTION
  int ii,i,k,j;
  dvar_vector preds_sexr(styr,endyr);
  dvar_matrix tmpo(1,2,styr,endyr);
  dvar_matrix tmpp(1,2,styr,endyr);
  dvar_vector obs_tmp(styr,endyr);
  dvariable ghl,ghl_number;
  dvariable hrate;
    
  Misc_output();
  tmpp1=0.0;
  tmpp2=0.0;
  tmpp3=0.0;
  tmpp4=0.0;
  report << Fout << " "<<fmort_pen<<endl;
  for (i=styr;i<=endyr;i++)
    if((totn_srv1(1,i)+totn_srv1(2,i))<0.01) {
      preds_sexr(i)=0.0;}
    else{
      preds_sexr(i)=totn_srv1(1,i)/(totn_srv1(1,i)+totn_srv1(2,i));
    }
  report << "Estimated numbers of immature new shell female crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) report <<  i<<" "<<natlength_inew(1,i) << endl;
  report << "Estimated numbers of immature old shell female crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) report <<  i<<" "<<natlength_iold(1,i) << endl;
  report << "Estimated numbers of mature new shell female crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) report <<  i<<" "<<natlength_mnew(1,i) << endl;
  report << "Estimated numbers of mature old shell female crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++)report <<  i<<" "<<natlength_mold(1,i) << endl;

  report << "Estimated numbers of immature new shell male crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) report << i<<" "<<natlength_inew(2,i) << endl;
  report << "Estimated numbers of immature old shell male crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) report << i<<" "<<natlength_iold(2,i) << endl;
  report << "Estimated numbers of mature new shell male crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) report << i<<" "<<natlength_mnew(2,i) << endl;
  report << "Estimated numbers of mature old shell male crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) report << i<<" "<<natlength_mold(2,i) << endl;
 
  report << "Observed numbers of immature new shell female crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(1,1,1,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  report << "Observed numbers of mature new shell female crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(2,1,1,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  report << "Observed numbers of mature old shell female crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(2,2,1,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  report << "Observed numbers of immature new shell male crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(1,1,2,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  report << "Observed numbers of immature old shell male crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(1,2,2,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  report << "Observed numbers of mature new shell male crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(2,1,2,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  report << "Observed numbers of mature old shell male crab by length: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(2,2,2,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  report << "Observed Survey Numbers by length females:  'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" " << obs_srv1_num(1,yrs_srv1_length(i)) << endl;
  report << "Observed Survey Numbers by length males: 'year', '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" " << obs_srv1_num(2,yrs_srv1_length(i))<< endl;
  
  report << "Predicted Survey Numbers by length females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" "  << pred_srv1(1,yrs_srv1_length(i)) << endl;
  report << "Predicted Survey Numbers by length males: 'year', '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) report<<yrs_srv1_length(i)<<" "  << pred_srv1(2,yrs_srv1_length(i)) << endl;
  report << "Predicted pop Numbers by length females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) report<<i<<" "<< natlength(1,i)<< endl;
  report << "Predicted pop Numbers by length males: 'year', '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) report<<i<<" "<< natlength(2,i)<< endl;
 
    //actual years for obs survey male are 1969,1970,1972-2009
   report<<"observed number of males greater than 101 mm: seq(1974,"<<endyr<<")"<<endl;
  report<<obs_lmales<<endl;
  report<<"observed biomass of males greater than 101 mm: seq(1974,"<<endyr<<")"<<endl;
  report<<obs_lmales_bio<<endl;
  report<<"pop estimate numbers of males >101: seq("<<styr<<","<<endyr<<")"<<endl;
  report<<legal_males<<endl;
  report<<"estimated population biomass of males > 101: seq("<<styr<<","<<endyr<<") "<<endl;
  report<<legal_males_bio<<endl;
  report<<"estimated survey numbers of males > 101: seq("<<styr<<","<<endyr<<") "<<endl;
  report<<legal_srv_males<<endl;
  report<<"estimated survey biomass of males > 101: seq("<<styr<<","<<endyr<<") "<<endl;
  report<<legal_srv_males_bio<<endl;
  report << "Observed survey biomass: seq(1974,"<<endyr<<")"<<endl;
  report << obs_srv1_biom(1974,endyr)<<endl;
  report << "predicted survey biomass: seq("<<styr<<","<<endyr<<")"<<endl;
  report << pred_srv1_bioms(1)+pred_srv1_bioms(2)<<endl;
  
  //survey numbers
  for(k=1;k<=2;k++)
   for(i=styr;i<=endyr;i++)
    {
     tmpo(k,i)=sum(obs_srv1_num(k,i));
     tmpp(k,i)=sum(pred_srv1(k,i));
    }
  report << "Observed survey numbers female: seq("<<styr<<","<<endyr<<")"<<endl;
  report << tmpo(1)<<endl;
  report << "Observed survey numbers male: seq("<<styr<<","<<endyr<<")"<<endl;
  report << tmpo(2)<<endl;
  report << "predicted survey numbers female: seq("<<styr<<","<<endyr<<")"<<endl;
  report << tmpp(1)<<endl;
  report << "predicted survey numbers male: seq("<<styr<<","<<endyr<<")"<<endl;
  report << tmpp(2)<<endl;
  report << "Observed survey female spawning biomass: seq("<<styr<<","<<endyr<<")"<<endl;
  report << obs_srv1_spbiom(1)<<endl;
  report << "Observed survey male spawning biomass: seq("<<styr<<","<<endyr<<")"<<endl;
  report << obs_srv1_spbiom(2)<<endl;
  report << "Observed survey female new spawning numbers: seq("<<styr<<","<<endyr<<")"<<endl;
  report << obs_srv1_spnum(1,1)<<endl;
  report << "Observed survey female old spawning numbers: seq("<<styr<<","<<endyr<<")"<<endl;
  report << obs_srv1_spnum(2,1)<<endl;
  report << "Observed survey male new spawning numbers: seq("<<styr<<","<<endyr<<")"<<endl;
  report << obs_srv1_spnum(1,2)<<endl;
  report << "Observed survey male old spawning numbers: seq("<<styr<<","<<endyr<<")"<<endl;
  report << obs_srv1_spnum(2,2)<<endl;
  report << "Observed survey female biomass: seq("<<styr<<","<<endyr<<")"<<endl;
  report << obs_srv1_bioms(1)<<endl;
  report << "Observed survey male biomass: seq("<<styr<<","<<endyr<<")"<<endl;
  report << obs_srv1_bioms(2)<<endl;
  report << "natural mortality immature females, males: 'FemM','MaleM'" << endl;
  report << M << endl;
  report << "natural mortality mature females, males: 'FemMm','MaleMm'" << endl;
  report << M_matn << endl;
  report << "natural mortality mature old shell females, males: 'FemMmo','MaleMmo'" << endl;
  report << M_mato << endl;
  report << "Predicted Biomass: seq("<<styr<<","<<endyr<<")" << endl;
  report << pred_bio << endl;
  report << "Predicted total population numbers: seq("<<styr<<","<<endyr<<") "<<endl;
  report <<popn<<endl;
  report << "Female Spawning Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio << endl;
  report << "Male Spawning Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << mspbio << endl;
  report << "Total Spawning Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio+mspbio << endl;
  report << "Female Spawning Biomass at fish time: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio_fishtime << endl;
  report << "Male Spawning Biomass at fish time: seq("<<styr<<","<<endyr<<") " << endl;
  report << mspbio_fishtime << endl;
  report << "Total Spawning Biomass at fish time: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio_fishtime+mspbio_fishtime << endl;
  report << "Mating time Female Spawning Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio_matetime << endl;
  report << "Mating time Male Spawning Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << mspbio_matetime << endl;
  report << "Mating time Male old shell Spawning Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << mspbio_old_matetime << endl;
  report << "Mating time female new shell Spawning Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio_new_matetime << endl;
  report << "Mating time Total Spawning Biomass : seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio_matetime+mspbio_matetime << endl;
  report << "Mating time effective Female Spawning Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << efspbio_matetime << endl;
  report << "Mating time effective Male Spawning Biomass(old shell only): seq("<<styr<<","<<endyr<<") " << endl;
  report << emspbio_matetime << endl;
  report << "Mating time Total effective Spawning Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << efspbio_matetime+emspbio_matetime << endl;
  report << "Mating time male Spawning numbers: seq("<<styr<<","<<endyr<<") " << endl;
  report << mspnum_matetime << endl;
  report << "Mating time Female Spawning numbers: seq("<<styr<<","<<endyr<<") " << endl;
  report << efspnum_matetime << endl;
  report << "Mating time Male Spawning numbers(old shell only): seq("<<styr<<","<<endyr<<") " << endl;
  report << emspnum_old_matetime << endl;
//  report << "ratio Mating time Female Spawning numbers to male old shell mature numbers : seq("<<styr<<","<<endyr<<") " << endl;
//  report << elem_div(efspnum_matetime,emspnum_old_matetime) << endl;
  report << "Mating time effective Female new shell Spawning biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report <<efspbio_new_matetime << endl;
  report << "Mating time Female new shell Spawning numbers: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspnum_new_matetime << endl;
//  report << "ratio Mating time Female new shell Spawning numbers to male old shell mature numbers : seq("<<styr<<","<<endyr<<") " << endl;
//            for(i=styr;i<=endyr;i++){if(emspnum_old_matetime(i)<0.001) emspnum_old_matetime(i)=1.0; }
//  report << elem_div(fspnum_new_matetime,emspnum_old_matetime) << endl;
  report << "Predicted Female survey Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << pred_srv1_bioms(1) << endl;
  report << "Predicted Male survey Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << pred_srv1_bioms(2)<< endl;
  report << "Predicted Female survey mature Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio_srv1 << endl;
  report << "Predicted Male survey mature Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << mspbio_srv1<< endl;
  report << "Predicted total survey mature Biomass: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio_srv1+mspbio_srv1<< endl;
  report << "Predicted Female survey new mature numbers: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio_srv1_num(1) << endl;
  report << "Predicted Female survey old mature numbers: seq("<<styr<<","<<endyr<<") " << endl;
  report << fspbio_srv1_num(2) << endl;
  report << "Predicted Male survey new mature numbers: seq("<<styr<<","<<endyr<<") " << endl;
  report << mspbio_srv1_num(1)<< endl;
  report << "Predicted Male survey old mature numbers: seq("<<styr<<","<<endyr<<") " << endl;
  report << mspbio_srv1_num(2)<< endl;

  report << "Observed Prop fishery ret new males:'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i<=nobs_fish; i++) report << yrs_fish(i) << " " << obs_p_fish_ret(1,i)<< endl;
  report << "Predicted length prop fishery ret new males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_fish; i++) 
   {
    ii=yrs_fish(i);  
    report <<  ii  <<  " "  <<  pred_p_fish_fit(1,ii)  << endl;
   }
  report << "Observed Prop fishery ret old males:'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i<=nobs_fish; i++) report << yrs_fish(i) << " " << obs_p_fish_ret(2,i)<< endl;
  report << "Predicted length prop fishery ret old males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_fish; i++)
   {
    ii=yrs_fish(i);  
    report <<  ii  <<  " "  <<  pred_p_fish_fit(2,ii)  << endl;
   }

  report << "Observed Prop fishery total new males:'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i<=nobs_fish_discm; i++) report << yrs_fish_discm(i) << " " << obs_p_fish_tot(1,i) << endl;
  report << "Predicted length prop fishery total new males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_fish_discm; i++)
   {
    ii=yrs_fish_discm(i);  
    report <<  ii  <<  " "  <<  pred_p_fish(1,ii)  << endl;
   }
  report << "Observed Prop fishery total old males:'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i<=nobs_fish_discm; i++) report << yrs_fish_discm(i) << " " << obs_p_fish_tot(2,i) << endl;
  report << "Predicted length prop fishery total old males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_fish_discm; i++)
   {
    ii=yrs_fish_discm(i);  
    report <<  ii  <<  " "  <<  pred_p_fish(2,ii)  << endl;
   }
  report << "Observed Prop fishery discard new males:'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i<=nobs_fish_discm; i++) report << yrs_fish_discm(i) << " " << obs_p_fish_discm(1,i) << endl;
  report << "Observed Prop fishery discard old males:'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for (i=1; i<=nobs_fish_discm; i++) report << yrs_fish_discm(i) << " " << obs_p_fish_discm(2,i)<< endl;

  report << "Observed length prop fishery discard all females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_fish_discf; i++) report <<  yrs_fish_discf(i)  <<  " "  <<  obs_p_fish_discf(i)  << endl;
  report << "Predicted length prop fishery discard all females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_fish_discf; i++)
   {
    ii=yrs_fish_discf(i);  
    report <<  ii  <<  " "  <<  pred_p_fish_discf(ii)  << endl;
   }
  report << "Observed length prop snow fishery females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_snowfish_discf; i++)
   {
    report <<  yrs_snowfish_discf(i)  <<  " "  <<  obs_p_snow(1,i)  << endl;
   }
  report << "Predicted length prop snow fishery females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_snowfish_discf; i++)
   {
    ii=yrs_snowfish_discf(i);  
    report <<  ii  <<  " "  <<  pred_p_snow(1,ii)  << endl;
   }
  report << "Observed length prop snow fishery males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_snowfish_discm; i++)
   {
    report <<  yrs_snowfish_discm(i)  <<  " "  <<  obs_p_snow(2,i)  << endl;
   }
  report << "Predicted length prop snow fishery males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_snowfish_discm; i++)
   {
    ii=yrs_snowfish_discm(i);  
    report <<  ii  <<  " "  <<  pred_p_snow(2,ii)  << endl;
   }
  report << "Observed length prop redk fishery females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_rkfish_discf; i++)
   {
    report <<  yrs_rkfish_discf(i)  <<  " "  <<  obs_p_rk(1,i)  << endl;
   }
  report << "Predicted length prop redk fishery females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_rkfish_discf; i++)
   {
    ii=yrs_rkfish_discf(i);  
    report <<  ii  <<  " "  <<  pred_p_rk(1,ii)  << endl;
   }
  report << "Observed length prop redk fishery males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_rkfish_discm; i++)
   {
    report <<  yrs_rkfish_discm(i)  <<  " "  <<  obs_p_rk(2,i)  << endl;
   }
  report << "Predicted length prop redk fishery males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_rkfish_discm; i++)
   {
    ii=yrs_rkfish_discm(i);  
    report <<  ii  <<  " "  <<  pred_p_rk(2,ii)  << endl;
   }

  report << "Predicted length prop trawl females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_trawl; i++)
   {
    ii=yrs_trawl(i);  
    report <<  ii  <<  " "  <<  pred_p_trawl(1,ii)  << endl;
   }
  report << "Observed length prop trawl females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_trawl; i++) report <<  yrs_trawl(i)  <<  " "  <<  obs_p_trawl(1,i)  << endl;
  report << "Predicted length prop trawl males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_trawl; i++)
   {
    ii=yrs_trawl(i);  
    report <<  ii  <<  " "  <<  pred_p_trawl(2,ii)  << endl;
   }
  report << "Observed length prop trawl males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_trawl; i++) report <<  yrs_trawl(i)  <<  " "  <<  obs_p_trawl(2,i)  << endl;

  report << "Observed Length Prop survey immature new females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    report << ii <<" " <<obs_p_srv1_len(1,1,1,i) << endl;
   }
  report << "Predicted length prop survey immature new females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    report << ii << " " << pred_p_srv1_len_new(1,1,ii) << endl;
   }
  report << "Observed Length Prop survey immature old females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    report << ii <<" " <<obs_p_srv1_len(1,2,1,i) << endl;
   }
  report << "Predicted length prop survey immature old females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    report << ii << " " << pred_p_srv1_len_old(1,1,ii) << endl;
   }
 
  report << "Observed Length Prop survey immature new males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++) report << yrs_srv1_length(i) <<" " <<obs_p_srv1_len(1,1,2,i) << endl;
  report << "Predicted length prop survey immature new males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    report << ii << " " << pred_p_srv1_len_new(1,2,ii) << endl;
   }
  report << "Observed Length Prop survey immature old males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++) report << yrs_srv1_length(i) <<" " <<obs_p_srv1_len(1,2,2,i) << endl;
  report << "Predicted length prop survey immature old males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
  {
   ii=yrs_srv1_length(i);  
   report << ii << " " << pred_p_srv1_len_old(1,2,ii) << endl;
  }
  report << "Observed Length Prop survey mature new females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    report << ii <<" " <<obs_p_srv1_len(2,1,1,i) << endl;
   }
  report << "Predicted length prop survey mature new females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    report << ii << " " << pred_p_srv1_len_new(2,1,ii) << endl;
   }
  report << "Observed Length Prop survey mature old females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    report << ii <<" " <<obs_p_srv1_len(2,2,1,i) << endl;
   }
  report << "Predicted length prop survey mature old females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    report << ii << " " << pred_p_srv1_len_old(2,1,ii) << endl;
   }
 
  report << "Observed Length Prop survey mature new males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++) report << yrs_srv1_length(i) <<" " <<obs_p_srv1_len(2,1,2,i) << endl;
  report << "Predicted length prop survey mature new males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    report << ii << " " << pred_p_srv1_len_new(2,2,ii) << endl;
   }
  report << "Observed Length Prop survey mature old males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++) report << yrs_srv1_length(i) <<" " <<obs_p_srv1_len(2,2,2,i) << endl;
  report << "Predicted length prop survey mature old males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    report << ii << " " << pred_p_srv1_len_old(2,2,ii) << endl;
   }
//for females don't have length data in first four years first year is 1974
     report << "Observed Length Prop survey all females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    report << ii <<" " <<obs_p_srv1_len(1,1,1,i)+obs_p_srv1_len(2,1,1,i)+obs_p_srv1_len(2,2,1,i)<< endl;
              tmpp4+=obs_p_srv1_len(1,1,1,i)+obs_p_srv1_len(2,1,1,i)+obs_p_srv1_len(2,2,1,i);
   }
  report << "Predicted length prop survey all females: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    report << ii << " " << pred_p_srv1_len_new(1,1,ii)+pred_p_srv1_len_new(2,1,ii)+pred_p_srv1_len_old(2,1,ii) << endl;
    tmpp1+=pred_p_srv1_len_new(1,1,ii)+pred_p_srv1_len_new(2,1,ii)+pred_p_srv1_len_old(2,1,ii);
   }
  report << "Observed Length Prop survey all males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    report << ii <<" " <<obs_p_srv1_len(1,1,2,i)+obs_p_srv1_len(1,2,2,i)+obs_p_srv1_len(2,1,2,i)+obs_p_srv1_len(2,2,2,i)<< endl;
         tmpp2+=obs_p_srv1_len(1,1,2,i)+obs_p_srv1_len(1,2,2,i)+obs_p_srv1_len(2,1,2,i)+obs_p_srv1_len(2,2,2,i);
   }
  report << "Predicted length prop survey all males: 'year','27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    report << ii << " " << pred_p_srv1_len_new(1,2,ii)+pred_p_srv1_len_new(2,2,ii)+pred_p_srv1_len_old(2,2,ii) << endl;
  tmpp3+=pred_p_srv1_len_new(1,2,ii)+pred_p_srv1_len_new(2,2,ii)+pred_p_srv1_len_old(2,2,ii);
   }
  report << "Sum of predicted prop survey all females: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
            report <<tmpp1<<endl;
  report << "Sum of predicted prop survey all males: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
            report <<tmpp3<<endl;
  report << "Sum of Observed prop survey all females: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
            report <<tmpp4<<endl;
  report << "Sum of Observed prop survey all males: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'" << endl;
            report <<tmpp2<<endl;

  report << "Predicted mean postmolt length females:  '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report << mean_length(1) << endl;
  report << "Predicted mean postmolt length males:  '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report << mean_length(2)<<endl; 
  report << "af1: 'females'" << endl;
  report << af1 << endl;
//  report << "af2: 'females'" << endl;
//  report << af2 << endl;
  report << "am1: 'males'" << endl;
  report << am1 << endl;
//  report << "am2: 'males'" << endl;
//  report << am2 << endl;
  report << "bf1: 'females'" << endl;
  report << bf1 << endl;
//  report << "bf2: 'females'" << endl;
//  report << bf2 << endl;
  report << "bm1: 'males'" << endl;
  report << bm1 << endl;
//  report << "bm2: 'males'" << endl;
//  report << bm2 << endl;
  report<<"Predicted probability of maturing females: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<<endl;
  report<<maturity_est(1)<<endl;
  report<<"Predicted probability of maturing males: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<<endl;
  report<<maturity_est(2)<<endl;
  report<<"molting probs female: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<<endl;
  report<<moltp(1)<<endl;
  report<<"molting probs male:'27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report<<moltp(2)<<endl;
  report <<"Molting probability mature males: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<moltp_mat(2)<<endl;
  report << "observed pot fishery cpue 1979 fishery to endyr fishery: seq(1970,"<<endyr<<")" << endl;
  report <<cpue(1969,endyr-1)<<endl;
  report << "predicted pot fishery cpue 1978 to endyr-1 survey: seq("<<styr+1<<","<<endyr<<")" << endl;
  report <<cpue_pred(styr,endyr-1)<<endl;
  report << "observed retained catch biomass: seq(1965,"<<endyr-1<<")" << endl;
  report << catch_ret(1965,endyr-1) << endl;
  report << "predicted retained catch biomass: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << pred_catch_ret(styr,endyr-1)<<endl;
  report << "predicted retained new catch biomass: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << (catch_male_ret_new*wtm)(styr,endyr-1)<<endl;
  report << "predicted retained old catch biomass: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << (catch_male_ret_old*wtm)(styr,endyr-1)<<endl;
  report << "observed retained+discard male catch biomass: seq(1992,"<<endyr-1<<")" << endl;
  report << obs_catchtot_biom(1992,endyr-1) << endl;
  report << "predicted retained+discard male catch biomass: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << pred_catch(styr,endyr-1) << endl;
  report << "predicted retained+discard new male catch biomass: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << (catch_lmale_new*wtm)(styr,endyr-1) << endl;
  report << "predicted retained+discard old male catch biomass: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << (catch_lmale_old*wtm)(styr,endyr-1) << endl;
  report << "observed discard male mortality biomass: seq(1992,"<<endyr-1<<")"<<endl;
  report << (obs_catchtot_biom(1992,endyr-1)-catch_ret(1992,endyr-1)) <<endl;
  report << "predicted discard male catch biomass: seq("<<styr<<","<<endyr-1<<")" << endl;
  report << pred_catch(styr,endyr-1) -pred_catch_ret(styr,endyr-1)<< endl;
  report << "observed female discard mortality biomass: seq(1992,"<<endyr-1<<")" << endl;
  report << obs_catchdf_biom(1992,endyr-1) << endl;
  report << "predicted female discard mortality biomass: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << pred_catch_disc(1)(styr,endyr-1) << endl;
  report << "observed male discard mortality biomass: seq(1992,"<<endyr-1<<")" << endl;
  report << obs_catchdm_biom(1992,endyr-1) << endl;
  report << "observed trawl catch biomass: seq("<<yrs_trawl_c(1)<<","<<yrs_trawl_c(nobs_trawl_c)<<")"<<endl;
  report << obs_catcht_biom(yrs_trawl_c)<<endl;
  report << "predicted trawl catch biomass: seq("<<styr<<","<<endyr<<")"<<endl;
  report <<pred_catch_trawl<<endl;
  report << "observed snow female discard mortality biomass: seq(1992,"<<endyr-1<<")" << endl;
   for (i=1; i<=nobs_discardc; i++)
    {
      report << catch_snowodisc(1)(i)<<" ";
     }
      report<< endl;
  report << "predicted snow female discard mortality biomass: seq("<<styr<<","<<endyr<<")" << endl;
  report << pred_catch_female_snowd << endl;
  report << "observed snow male discard mortality biomass: seq(1992,"<<endyr-1<<")" << endl;
   for (i=1; i<=nobs_discardc; i++)
    {
      report << catch_snowodisc(2)(i) <<" ";
     }
      report<< endl;
  report << "predicted snow male discard mortality biomass: seq("<<styr<<","<<endyr<<")" << endl;
  report << pred_catch_snowd << endl;  
  report << "observed redk female discard mortality biomass: seq(1992,"<<endyr-1<<")" << endl;
   for (i=1; i<=nobs_discardc; i++)
    {
      report << catch_rkodisc(1)(i) <<" ";
     }
  report << endl;
  report << "predicted redk female discard mortality biomass: seq("<<styr<<","<<endyr<<")" << endl;
  report << pred_catch_female_rkd << endl;
  report << "observed redk male discard mortality biomass: seq(1992,"<<endyr-1<<")" << endl;
   for (i=1; i<=nobs_discardc; i++)
    {
      report << catch_rkodisc(2)(i) <<" ";
     }
  report << endl;
  report << "predicted redk male discard mortality biomass: seq("<<styr<<","<<endyr<<")" << endl;
  report << pred_catch_rkd << endl;
  report << "predicted total male catch biomass: seq("<<styr+1<<","<<endyr<<")" << endl;
  report<<pred_catch(styr,endyr-1)+pred_catch_rkd(styr,endyr-1)+pred_catch_snowd(styr,endyr-1)+pred_catch_trawl(styr,endyr-1)/2.0<<endl;
  report << "predicted total female catch biomass: seq("<<styr+1<<","<<endyr<<")" << endl;
  report<<pred_catch_disc(1)(styr,endyr-1)+pred_catch_female_rkd(styr,endyr-1)+pred_catch_female_snowd(styr,endyr-1)+pred_catch_trawl(styr,endyr-1)/2.0<<endl;
//  report<<"Estimated total catch div. by male spawing biomass at fishtime: seq("<<styr+1<<","<<endyr<<")"<<endl;
//  report<<elem_div(pred_catch(styr,endyr-1)+pred_catch_rkd(styr,endyr-1)+pred_catch_snowd(styr,endyr-1)+pred_catch_trawl(styr,endyr-1)/2.0,mspbio_fishtime(styr,endyr-1))<<endl;
//  report << "estimated retained catch div. by male spawning biomass at fishtime: seq("<<styr+1<<","<<endyr<<")" << endl;
//  report <<elem_div(pred_catch_ret,mspbio_fishtime)(styr,endyr-1) << endl;
//  report << "estimated total catch div. by male spawning biomass at fishtime: seq("<<styr+1<<","<<endyr<<")" << endl;
//  report <<elem_div(pred_catch,mspbio_fishtime)(styr,endyr-1) << endl;
//  report << "estimated total catch of legal males by legal males at fishtime: seq("<<styr+1<<","<<endyr<<")" << endl;
//   report <<elem_div(pred_catch_gt101(styr,endyr-1),bio_males_gt101(styr,endyr-1)) << endl;
//  report << "estimated total catch numbers of males >101 div. by males numbers >101 at fishtime: seq(1970,"<<endyr<<")" << endl;
//  report <<elem_div(pred_catch_no_gt101(styr,endyr-1),num_males_gt101(styr,endyr-1)) << endl;
//  report << "estimated total catch numbers of males >101 div. by survey estimate males numbers >101 at fishtime: seq("<<styr+1<<","<<endyr<<")" << endl;
//  for(i=styr;i<endyr;i++) obs_tmp(i) = obs_lmales(i-(styr-1));
//  report <<elem_div(pred_catch_no_gt101(styr,endyr-1),obs_tmp(styr,endyr-1)*mfexp(-M_matn(2)*(7/12))) << endl;
//  for(i=styr;i<endyr;i++) obs_tmp(i) = obs_lmales_bio(i-(styr-1));
//  report << "estimated total catch biomass of males >101 div. by survey estimate male biomass >101 at fishtime: seq("<<styr+1<<","<<endyr<<")" << endl;
//  report <<elem_div(pred_catch_gt101(styr,endyr-1),obs_tmp(styr,endyr-1)*mfexp(-M_matn(2)*(7/12)) ) << endl;
//  report << "estimated total catch biomass div. by survey estimate male mature biomass at fishtime: seq("<<styr+1<<","<<endyr<<")" << endl;
//  report <<elem_div(pred_catch(styr,endyr-1),((obs_srv1_spbiom(2))(styr,endyr-1))*mfexp(-M_matn(2)*(7/12))) << endl;
  report << "estimated annual total directed fishing mortality: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << fmort(styr,endyr-1) << endl;
  report << "estimated annual snow fishing mortality: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << fmortd_snow(styr,endyr-1) << endl;
  report << "estimated annual red king fishing mortality: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << fmortd_rk(styr,endyr-1) << endl;
  report << "estimated annual total fishing mortality: seq("<<styr+1<<","<<endyr<<")" << endl;
  for(i=styr;i<=(endyr-1);i++) report << F(1,i)(nlenm)+ Fdisct(2,i)(nlenm)+Fdisc_snow(2,i)(nlenm)+Fdisc_rk(2,i)(nlenm) <<" "; report<< endl;
  report <<"retained F: seq("<<styr+1<<","<<endyr<<")" << endl;
  for(i=styr;i<=(endyr-1);i++) report <<F_ret(1,i)(nlenm)<<" "; report<<endl;
  report <<"ghl: seq(1979,"<<endyr<<")" << endl;
  report <<catch_ghl/2.2<<endl;
  report << "estimated annual fishing mortality females pot: seq("<<styr+1<<","<<endyr<<")" << endl;
  for(i=styr;i<=(endyr-1);i++) report << Fdiscf(i)(nlenm) <<" "; report<<endl;
  report << "estimated annual fishing mortality trawl bycatch: seq("<<styr+1<<","<<endyr<<")" << endl;
  report << fmortt(styr,endyr-1) <<endl;
//recruits in the model are 1978 to 2004, the 1978 recruits are those that enter the population
//in spring of 1979, before the 1979 survey - since using the survey as the start of the year
// in the model spring 1979 is stil 1978.  the last recruits are 2003 that come in spring 2004
  report << "estimated number of recruitments female: seq("<<styr+1<<","<<endyr<<")" << endl;
  for(i=styr; i<=1973; i++) report << mfexp(mean_log_rec1_early+rec_devf_early(i))<<" ";
  for(i=1974; i<=(endyr-1); i++) report << mfexp(mean_log_rec(1)+rec_dev(1,i))<<" ";
  report <<endl<< "estimated number of recruitments male: seq("<<styr+1<<","<<endyr<<")" << endl;
  for(i=styr; i<=1973; i++) report << mfexp(mean_log_rec1_early+rec_devf_early(i))<<" ";
  for(i=1974; i<=(endyr-1); i++) report << mfexp(mean_log_rec(2)+rec_dev(2,i))<<" ";
  for(i=1;i<=median_rec_yrs;i++)report<<2*median_rec<<" "; report<<endl;
  
  report<<"distribution of recruits to length bins: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<<endl;
  report<<rec_len<<endl;
//  report<<"fishery total selectivity new shell 50% parameter: seq("<<1981+1<<","<<endyr<<")"<<endl;
//  report <<fish_sel50_mn(1981,endyr-1)<<endl;
//  report <<"fishery total selectivity old shell 50% parameter: seq("<<styr+1<<","<<endyr<<")"<<endl;
//  report <<mfexp(log_avg_sel50_mo+log_sel50_dev_mo)(styr,endyr-1)<<endl;
  report << "selectivity fishery total new males styr to 1991: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report << sel(1,1990) << endl;
//  report << "selectivity fishery total new males 1992 to 1996: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
//  report << sel(1,1991) << endl;
  report << "selectivity fishery total new males 1991 to 1996: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
   for(i=1991;i<=1996; i++) 
   { report << sel(1,i) << endl;}
  report << "selectivity fishery total new males 2005 to present: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
   for(i=2005;i<=endyr; i++) 
   { report << sel(1,i) << endl;}
//  report << "selectivity fishery total old males: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
//  report << sel(2) << endl;
  report << "selectivity fishery ret new males styr to 1991: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;  
  report << sel_fit(1,1990) << endl;
//  report << "selectivity fishery ret new males 1992 to 1996: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;  
//  report << sel_fit(1,1991) << endl;
  report << "selectivity fishery ret new males 1991 to 1996: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;  
  for(i=1991;i<=1996; i++) 
   { for(j=1; j<=nlenm; j++){report <<sel_fit(1,i,j)<<" ";
   }
      report<<endl;
   }
  report << "selectivity fishery ret new males 2005 to present: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;  
  for(i=2005;i<=endyr; i++) 
   { for(j=1; j<=nlenm; j++){report <<sel_fit(1,i,j)<<" ";
   }
      report<<endl;
   }
//  report << "selectivity fishery ret new males: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
//  report << sel_fit(1) << endl;
//  report << "selectivity fishery ret old males: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
//  report << sel_fit(2) << endl;
//  report <<"retention curve males new: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
//  report <<sel_ret<<endl;
  report <<"retention curve males new: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<sel_ret(1,endyr-1)<<endl;
  report <<"retention curve males old: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<sel_ret(2,endyr-1)<<endl;
  report << "selectivity discard females: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<sel_discf<<endl;
  report << "selectivity trawl females:'27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<sel_trawl(1,1)<<endl;
  report <<sel_trawl(2,1)<<endl;
  report <<sel_trawl(3,1)<<endl;
  report << "selectivity trawl males:'27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<sel_trawl(1,2)<<endl;
  report <<sel_trawl(2,2)<<endl;
  report <<sel_trawl(3,2)<<endl;
  report << "selectivity snow females:'27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<sel_disc_snow(1,1)<<endl;
  report <<sel_disc_snow(2,1)<<endl;
  report <<sel_disc_snow(3,1)<<endl;
  report << "selectivity snow males:'27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<sel_disc_snow(1,2)<<endl;
  report <<sel_disc_snow(2,2)<<endl;
  report <<sel_disc_snow(3,2)<<endl;
  report << "selectivity redk females:'27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<sel_disc_rkc(1,1)<<endl;
  report <<sel_disc_rkc(2,1)<<endl;
  report <<sel_disc_rkc(3,1)<<endl;
  report << "selectivity redk males:'27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report <<sel_disc_rkc(1,2)<<endl;
  report <<sel_disc_rkc(2,2)<<endl;
  report <<sel_disc_rkc(3,2)<<endl;
//  report << "selectivity survey females 1969 1973: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
//  report << sel_srv1(1) << endl;
//  report << "selectivity survey males 1969 1973: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
//  report << sel_srv1(2) << endl;
  report << "selectivity survey females 1974 to 1981: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report << "selectivity survey males 1974 to 1981: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report << "selectivity survey females 1982 to 1987: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report << "selectivity survey males 1982 to 1987: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report << "selectivity survey females 1988 to endyr: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report << "selectivity survey males 1988 to endyr: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report << "numbers of mature females by age and length: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++)
   { report << natlength_mnew(1,i)<<endl; report << natlength_mold_age(1,i)<<endl; }
  
  report << "numbers of mature males by age and length: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  for(i=styr;i<=endyr;i++) 
   { report << natlength_mnew(2,i)<<endl; report << natlength_mold_age(2,i)<<endl; }
  report << "pred_sexr population: seq("<<styr<<","<<endyr<<")" << endl;
  report << predpop_sexr <<endl;
  report << "pred_sexr survey: seq("<<styr<<","<<endyr<<")" << endl;
  report << preds_sexr <<endl;
  report <<"likelihood: 'penal_rec','penal_sexr','change_sel_like','len_like_ret','len_like_tot','len_like_fem','len_like_surv','len_like_trawl', 'fpen',  'catch_like tot','catch ret', 'catch fem','catch trawl', 'surv_like','surv_like_nowt','largemale_like','initnum_penal','initsmo_penal','total likelihood'"<<endl;
  report <<penal_rec<<"  "<<penal_sexr<<"  "<<sel_50m_penal<<"  "
           <<like_wght(1)*len_like(1)<<" "<<like_wght(2)*len_like(2)<<" "<<like_wght(3)*len_like(3)<<" "<<like_wght(4)*len_like(4)<<" "<<like_wght(7)*len_like(5)<<" "<< " "
           <<fpen<<" "<<wght_total_catch*catch_like1<<" "<<like_wght(6)*catch_like2<<" "<<wght_female_potcatch*catch_likef<<" "<<like_wght(6)*0.01*catch_liket<<" "<<surv_like<<" "<<surv_like_nowt<<" "<<wt_lmlike*largemale_like<<" "<<initnum_penal<<" "<<initsmo_penal<<" "<<f<<endl;
  report <<"offset for survey lengths"<<endl;
  report <<offset(4)<<endl;
  report <<"survey length likelihoods: 'immature new female','immature new male','immature old female','immature old male','mature new female','mature new male','mature old female','mature old male'"<<endl;
  report <<len_like_srv<<endl;
  report <<"likelihood weights: 'retained length','total catch length','female catch','survey length','survey biomass','catch biomass','trawl length'"<<endl;
  report <<like_wght<<"  "<<like_wght_mbio<<endl;
  report <<"likelihood weights:  'rec devs','sex ratio','fishery 50%','fmort phase 1','fmort phase>1','fmort devs'"<<endl;
  report <<like_wght_rec<<"  "<<like_wght_sexr<<"  "<<like_wght_sel50<<"  "<<like_wght_fph1<<"  "<<like_wght_fph2<<"  "<<like_wght_fdev<<endl;
  report <<"likes bayesian: 'like_mmat','af_penal','bf_penal','am_penal','bm_penal'"<<endl;
  report <<af_penal<<" "<<bf_penal<<" "<<am_penal<<" "<<bm_penal<<endl;
  report<<"length - length transition matrix Females:'27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report<<len_len(1)<<endl;  
  report<<"length - length transition matrix Males:'27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;
  report<<len_len(2)<<endl;  
  report<<"effective N survey lengths immature new shell female "<<endl;
  report<<effn_srv1(1,1,1)<<endl;
  report<<" effective N survey lengths mature new shell female "<<endl;
  report<<effn_srv1(2,1,1)<<endl;
  report<<" effective N survey lengths mature old shell female "<<endl;
  report<<effn_srv1(2,2,1)<<endl;
  report<<" effective N survey lengths immature new shell male "<<endl;
  report<<effn_srv1(1,1,2)<<endl;
  report<<" effective N survey lengths immature old shell male "<<endl;
  report<<effn_srv1(1,2,2)<<endl;
  report<<" effective N survey lengths mature new shell male "<<endl;
  report<<effn_srv1(2,1,2)<<endl;
  report<<" effective N survey lengths mature old shell male "<<endl;
  report<<effn_srv1(2,2,2)<<endl;
  report<<" effective N retained lengths new, old shell "<<endl;
  report<<effn_fish_ret<<endl;
  report<<" effective N total lengths new, old shell"<<endl;
  report<<effn_fish_tot<<endl;
  report<<"male new shell total pot fishery exploitation rates"<<endl;
  report<<1-mfexp(-1.0*F(1))<<endl;
  report<<"male old shell total pot fishery exploitation rates"<<endl;
  report<<1-mfexp(-1.0*F(2))<<endl;
  report<<"numbers new shell males at time of pop fishery"<<endl;
  report<<natl_new_fishtime(2)<<endl;
  report<<"numbers old shell males at time of pop fishery"<<endl;
  report<<natl_old_fishtime(2)<<endl;
  report<<"total catch in numbers new shell males"<<endl;
  report<<catch_lmale_new<<endl;
  report<<"total catch in numbers old shell males"<<endl;
  report<<catch_lmale_old<<endl;
  report<<"retained catch in numbers new shell males"<<endl;
  report<<catch_male_ret_new<<endl;
  report<<"retained catch in numbers old shell males"<<endl;
  report<<catch_male_ret_old<<endl;
  report<<"observed retained catch new shell males"<<endl;
  for (i=1; i<=nobs_fish; i++) 
   report<<yrs_fish(i)<<" "<<obs_p_fish_ret(1,i)*catch_numbers(yrs_fish(i))<<endl;
  report<<"observed retained catch old shell males"<<endl;
  for (i=1; i<=nobs_fish; i++) 
   report<<yrs_fish(i)<<" "<<obs_p_fish_ret(2,i)*catch_numbers(yrs_fish(i))<<endl;
  
  //compute GHL 
  for(i=2000;i<=endyr;i++)
   {

    hrate=0.1+(((mspbio_srv1(i)+fspbio_srv1(i))*2.2-230.4)*(0.125/691.2));
    if((mspbio_srv1(i)+fspbio_srv1(i))<=(230.4/2.2)) hrate=0.0;
    if((mspbio_srv1(i)+fspbio_srv1(i))>(921.6/2.2)) hrate=0.225;
    ghl=hrate*2.2*mspbio_srv1(i);
    // get numbers by dividing by average weight of crabs greater than 102 from 2003 survey
    ghl_number=ghl/1.27;
    // cap of 58% of exploitable males = new shell>101 + 25% of old shell>101
    if((1000.*ghl_number)> (0.58*(legal_srv_males_n(i)+(0.25*legal_srv_males_o(i)))))
     {
      ghl_number = (0.58*(legal_srv_males_n(i)+(0.25*legal_srv_males_o(i))))/1000.;
      ghl = ghl_number*1.27;
     }

    report <<"year, harvest rate, GHL in 1000 tons then 1000's of crabs: 'year','hrate','ghl','ghl_nos'";
    report <<i<<"  "<<hrate<<"  "<<ghl<<"  "<<ghl_number<<endl;
    report <<" estimated survey mature female biomass "<<fspbio_srv1(i)<<endl;
    report <<" estimated survey mature male biomass "<<mspbio_srv1(i)<<endl;
    report <<"number of estimated survey new males > 101mm "<<legal_srv_males_n(i)<<endl;
    report <<"number of estimated survey old males > 101mm "<<legal_srv_males_o(i)<<endl;

   }
  // stuff for input to projection model
  report<<"#number of length bins"<<endl;
  report<<nlenm<<endl;
  report<<"#Nat mort immature female/male"<<endl;
  report<<M<<endl;
  report<<"#nat mort mature new shell female/male"<<endl;
  report<<M_matn<<endl;
  report<<"#nat mort mature old shell female/male"<<endl;
  report<<M_mato<<endl;
  report<<"#constant recruitment"<<endl;
  report<<"1000000"<<endl;
  report<<"#average of last 4 years sel total male new old shell"<<endl;
  report<<(sel(1,endyr-4)+sel(1,endyr-3)+sel(1,endyr-2)+sel(1,endyr-1))/4.0<<endl;
  report<<(sel(1,endyr-4)+sel(2,endyr-3)+sel(2,endyr-2)+sel(2,endyr-1))/4.0<<endl;
  report<<"#average of last 4 years sel retained curve male new old shell"<<endl;
  report<<(sel_fit(1,endyr-3)+sel_fit(1,endyr-2)+sel_fit(1,endyr-1))/3.0<<endl;
  report<<(sel_fit(2,endyr-3)+sel_fit(2,endyr-2)+sel_fit(2,endyr-1))/3.0<<endl;
  report<<"#trawl selectivity female male"<<endl;
  report<<sel_trawl(3)<<endl;
  report<<"#female pot discard selectivity"<<endl;
  report<<sel_discf<<endl;
  report <<"#selectivity snow females"<< endl;
  report <<sel_disc_snow(3,1)<<endl;
  report << "#selectivity snow males"<< endl;
  report <<sel_disc_snow(3,2)<<endl;
  report <<"#selectivity redk females"<< endl;
  report <<sel_disc_rkc(3,1)<<endl;
  report <<"#selectivity redk males"<< endl;
  report <<sel_disc_rkc(3,2)<<endl; 
  report<<"#maturity curve new shell female male"<<endl;
  report<<maturity_est(1)<<endl;
  report<<maturity_est(2)<<endl;
  report<<"#maturity curve old shell female male"<<endl;
  report<<maturity_old_average<<endl;
  report<<"#molting probability immature female male"<<endl;
  report<<moltp<<endl;
  report<<"#molting probability mature female male"<<endl;
  report<<moltp_mat<<endl;
  report<<"#prop recruits to new shell"<<endl;
  report<<proprecn<<endl;
  report<<"#distribution of recruits to length bins"<<endl;
  report<<rec_len<<endl;
  report<<"#time of catch in fraction of year from survey - 7 months"<<endl;
  report<<catch_midpt(endyr)<<endl;
  report<<"#number at length new shell females males at time of fishery endyr from model"<<endl;
  report<<natl_new_fishtime(1,endyr)<<endl;
  report<<natl_new_fishtime(2,endyr)<<endl;
  report<<"#number at length old shell females males at time of fishery endyr from model"<<endl;
  report<<natl_old_fishtime(1,endyr)<<endl;
  report<<natl_old_fishtime(2,endyr)<<endl;
  report<<"#last year male spawning biomass"<<endl;
  report<<mspbio(endyr)<<endl;
  report<<"#last year female spawning biomass"<<endl;
  report<<fspbio(endyr)<<endl;
  report<<"#last year male spawning biomass at matingtime"<<endl;
  report<<mspbio_matetime(endyr)<<endl;
  report<<"#last year female spawning biomass at matingtime"<<endl;
  report<<fspbio_matetime(endyr)<<endl;
  report<<"#numbers at length immature new shell female male last year"<<endl;
  report<<natlength_inew(1,endyr)<<endl;
  report<<natlength_inew(2,endyr)<<endl;
  report<<"#numbers at length immature old shell female male last year"<<endl;
  report<<natlength_iold(1,endyr)<<endl;
  report<<natlength_iold(2,endyr)<<endl;
  report<<"#numbers at length mature new shell female male last year"<<endl;
  report<<natlength_mnew(1,endyr)<<endl;
  report<<natlength_mnew(2,endyr)<<endl;
  report<<"#numbers at length mature old shell female male last year"<<endl;
  report<<natlength_mold(1,endyr)<<endl;
  report<<natlength_mold(2,endyr)<<endl;
  report<<"#weight at length female juvenile"<<endl;
  report<<wtf(1)<<endl;
  report<<"#weight at length female mature"<<endl;
  report<<wtf(2)<<endl;
  report<<"#weight at length male"<<endl;
  report<<wtm<<endl;
  report<<"#length-length transition matrix"<<endl;
  report<<len_len<<endl;
  report<<"#female discard pot fishing F"<<endl;
  report<<mean(fmortdf(endyr-5,endyr-1))<<endl;
  report<<"#trawl fishing F female male average last 5 yrs"<<endl;
  report<<mean(fmortt(endyr-5,endyr-1))<<endl;
  report<<"#snow fishing F female male average last 5 yrs"<<endl;
  report<<mean(fmortd_snow(endyr-5,endyr-1))<<endl;
  report<<"#red king fishing F female male average last 5 yrs"<<endl;
  report<<mean(fmortd_rk(endyr-5,endyr-1))<<endl;
  report<<"#number of recruits from the model styr to endyr-1"<<endl;
  report<<endyr-styr<<endl;
  report<<"#number of recruits for avg to estimate B35%"<<endl;
  report<<endyr-1960<<endl;
  report <<"#recruitments female, male start year to endyr-1 from model" << endl;
  for(i=styr; i<=1973; i++) report << mfexp(mean_log_rec1_early+rec_devf_early(i))<<" ";
  for(i=1974; i<=(endyr-1); i++) report << mfexp(mean_log_rec(1)+rec_dev(1,i))<<" ";
  report <<endl<< "#recruitments male, male start 1960 to endyr-1 from model" << endl;
  for(i=styr; i<=1973; i++) report << mfexp(mean_log_rec1_early+rec_devf_early(i))<<" ";
  for(i=1974; i<endyr; i++) report << mfexp(mean_log_rec(2)+rec_dev(2,i))<<" "; report<<endl;
  
  report<<"#male spawning biomass at matetime for endyr-5 to endyr-1 for spawner recruit curve to estimate recruitments"<<endl;
  report<<mspbio_matetime(endyr-5,endyr-1)<<endl;
  report<<"#male spawning biomass at matetime for str year to endyr-1 for spawner recruit curve to estimate recruitments"<<endl;
  report<<mspbio_matetime(styr,endyr-1)<<endl;
  report <<"#selectivity survey males 1989 to endyr: '27.5','32.5','37.5','42.5','47.5','52.5','57.5','62.5','67.5','72.5','77.5','82.5','87.5','92.5','97.5','102.5','107.5','112.5','117.5','122.5','127.5','132.5','137.5','142.5','147.5','152.5','157.5','162.5','167.5','172.5','177.5','182.5'"<< endl;


  if (last_phase())
  {
    R_out << "$Estimated.numbers.of.immature.new.shell.female.crab.by.length"<< endl;
    for(i=styr;i<=endyr;i++) R_out <<  i<<" "<<natlength_inew(1,i) << endl;
    R_out << "$Estimated.numbers.of.immature.old.shell.female.crab.by.length"<< endl;
    for(i=styr;i<=endyr;i++) R_out <<  i<<" "<<natlength_iold(1,i) << endl;
    R_out << "$Estimated.numbers.of.mature.new.shell.female.crab.by.length"<< endl;
    for(i=styr;i<=endyr;i++) R_out <<  i<<" "<<natlength_mnew(1,i) << endl;
    R_out << "$Estimated.numbers.of.mature.old.shell.female.crab.by.length"<< endl;
    for(i=styr;i<=endyr;i++)R_out <<  i<<" "<<natlength_mold(1,i) << endl;
  
    R_out << "$Estimated.numbers.of.immature.new.shell.male.crab.by.length"<< endl;
    for(i=styr;i<=endyr;i++) R_out << i<<" "<<natlength_inew(2,i) << endl;
    R_out << "$Estimated.numbers.of.immature.old.shell.male.crab.by.length"<< endl;
    for(i=styr;i<=endyr;i++) R_out << i<<" "<<natlength_iold(2,i) << endl;
    R_out << "$Estimated.numbers.of.mature.new.shell.male.crab.by.length"<< endl;
    for(i=styr;i<=endyr;i++) R_out << i<<" "<<natlength_mnew(2,i) << endl;
    R_out << "$Estimated.numbers.of.mature.old.shell.male.crab.by.length"<< endl;
    for(i=styr;i<=endyr;i++) R_out << i<<" "<<natlength_mold(2,i) << endl;
   
    R_out << "$Observed.numbers.of.immature.new.shell.female.crab.by.length"<< endl;
    for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(1,1,1,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
    R_out << "$Observed.numbers.of.mature.new.shell.female.crab.by.length"<< endl;
    for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(2,1,1,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
    R_out << "$Observed.numbers.of.mature.old.shell.female.crab.by.length"<< endl;
    for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(2,2,1,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
    R_out << "$Observed.numbers.of.immature.new.shell.male.crab.by.length"<< endl;
    for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(1,1,2,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  R_out << "$Observed.numbers.of.immature.old.shell.male.crab.by.length"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(1,2,2,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  R_out << "$Observed.numbers.of.mature.new.shell.male.crab.by.length"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(2,1,2,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  R_out << "$Observed.numbers.of.mature.old.shell.male.crab.by.length"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" "<<obs_p_srv1_len(2,2,2,i)*obs_srv1t(yrs_srv1_length(i))<<endl;
  R_out << "$Observed.Survey.Numbers.by.length.females"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" " << obs_srv1_num(1,yrs_srv1_length(i)) << endl;
  R_out << "$Observed.Survey.Numbers.by.length.males"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" " << obs_srv1_num(2,yrs_srv1_length(i))<< endl;
  
  R_out << "$Predicted.Survey.Numbers.by.length.females"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" "  << pred_srv1(1,yrs_srv1_length(i)) << endl;
  R_out << "$Predicted.Survey.Numbers.by.length.males"<< endl;
  for (i=1; i <= nobs_srv1_length; i++) R_out<<yrs_srv1_length(i)<<" "  << pred_srv1(2,yrs_srv1_length(i)) << endl;
  R_out << "$Predicted.pop.Numbers.by.length.females"<< endl;
  for(i=styr;i<=endyr;i++) R_out<<i<<" "<< natlength(1,i)<< endl;
  R_out << "$Predicted.pop.Numbers.by.length.males"<< endl;
  for(i=styr;i<=endyr;i++) R_out<<i<<" "<< natlength(2,i)<< endl;
 
    //actual years for obs survey male are 1969,1970,1972-2009
   R_out<<"$observed.number.of.males.greater.than.101.mm"<<endl;
  R_out<<obs_lmales<<endl;
  R_out<<"$observed.biomass.of.males.greater.than.101.mm"<<endl;
  R_out<<obs_lmales_bio<<endl;
  R_out<<"$pop.estimate.numbers.of.males.101"<<endl;
  R_out<<legal_males<<endl;
  R_out<<"$estimated.population.biomass.of.males.101"<<endl;
  R_out<<legal_males_bio<<endl;
  R_out<<"$estimated.survey.numbers.of.males.101"<<endl;
  R_out<<legal_srv_males<<endl;
  R_out<<"$estimated.survey.biomass.of.males.101"<<endl;
  R_out<<legal_srv_males_bio<<endl;
  R_out << "$Observed.survey.biomass"<<endl;
  R_out << obs_srv1_biom(1969,endyr)<<endl;
  R_out << "$predicted.survey.biomass"<<endl;
  R_out << pred_srv1_bioms(1)+pred_srv1_bioms(2)<<endl;
  
  //survey numbers
  for(k=1;k<=2;k++)
   for(i=styr;i<=endyr;i++)
    {
     tmpo(k,i)=sum(obs_srv1_num(k,i));
     tmpp(k,i)=sum(pred_srv1(k,i));
    }
  R_out << "$Observed.survey.numbers.female"<<endl;
  R_out << tmpo(1)<<endl;
  R_out << "$Observed.survey.numbers.male"<<endl;
  R_out << tmpo(2)<<endl;
  R_out << "$predicted.survey.numbers.female"<<endl;
  R_out << tmpp(1)<<endl;
  R_out << "$predicted.survey.numbers.male"<<endl;
  R_out << tmpp(2)<<endl;
  R_out << "$Observed.survey.female.spawning.biomass"<<endl;
  R_out << obs_srv1_spbiom(1)<<endl;
  R_out << "$Observed.survey.male.spawning.biomass"<<endl;
  R_out << obs_srv1_spbiom(2)<<endl;
  R_out << "$Observed.survey.female.new.spawning.numbers"<<endl;
  R_out << obs_srv1_spnum(1,1)<<endl;
  R_out << "$Observed.survey.female.old.spawning.numbers"<<endl;
  R_out << obs_srv1_spnum(2,1)<<endl;
  R_out << "$Observed.survey.male.new.spawning.numbers"<<endl;
  R_out << obs_srv1_spnum(1,2)<<endl;
  R_out << "$Observed.survey.male.old.spawning.numbers"<<endl;
  R_out << obs_srv1_spnum(2,2)<<endl;
  R_out << "$Observed.survey.female.biomass"<<endl;
  R_out << obs_srv1_bioms(1)<<endl;
  R_out << "$Observed.survey.male.biomass"<<endl;
  R_out << obs_srv1_bioms(2)<<endl;
  R_out << "$natural.mortality.immature.females.males" << endl;
  R_out << M << endl;
  R_out << "$natural.mortality.mature.females.males" << endl;
  R_out << M_matn << endl;
  R_out << "$natural.mortality.mature.old.shell.females.males" << endl;
  R_out << M_mato << endl;
  R_out << "$Predicted.Biomass" << endl;
  R_out << pred_bio << endl;
  R_out << "$Predicted.total.population.numbers"<<endl;
  R_out <<popn<<endl;
  R_out << "$Female.Spawning.Biomass" << endl;
  R_out << fspbio << endl;
  R_out << "$Male.Spawning.Biomass" << endl;
  R_out << mspbio << endl;
  R_out << "$Total.Spawning.Biomass" << endl;
  R_out << fspbio+mspbio << endl;
  R_out << "$Female.Spawning.Biomass.at.fish.time" << endl;
  R_out << fspbio_fishtime << endl;
  R_out << "$Male.Spawning.Biomass.at.fish.time" << endl;
  R_out << mspbio_fishtime << endl;
  R_out << "$Total.Spawning.Biomass.at.fish.time" << endl;
  R_out << fspbio_fishtime+mspbio_fishtime << endl;
  R_out << "$Mating.time.Female.Spawning.Biomass" << endl;
  R_out << fspbio_matetime << endl;
  R_out << "$Mating.time.Male.Spawning.Biomass" << endl;
  R_out << mspbio_matetime << endl;
  R_out << "$Mating.time.Male.old.shell.Spawning.Biomasss" << endl;
  R_out << mspbio_old_matetime << endl;
  R_out << "$Mating.time.female.new.shell.Spawning.Biomass" << endl;
  R_out << fspbio_new_matetime << endl;
  R_out << "$Mating.time.Total.Spawning.Biomass" << endl;
  R_out << fspbio_matetime+mspbio_matetime << endl;
  R_out << "$Mating.time.effective.Female.Spawning.Biomass" << endl;
  R_out << efspbio_matetime << endl;
  R_out << "$Mating.time.effective.Male.Spawning.Biomass.old.shell.only" << endl;
  R_out << emspbio_matetime << endl;
  R_out << "$Mating.time.Total.effective.Spawning.Biomass" << endl;
  R_out << efspbio_matetime+emspbio_matetime << endl;
  R_out << "$Mating.time.male.Spawning.numbers" << endl;
  R_out << mspnum_matetime << endl;
  R_out << "$Mating.time.Female.Spawning.numbers" << endl;
  R_out << efspnum_matetime << endl;
  R_out << "$Mating.time.Male.Spawning.numbers.old.shell.only" << endl;
  R_out << emspnum_old_matetime << endl;
//  R_out << "$ratio.Mating.time.Female.Spawning.numbers.to.male.old.shell.mature.numbers" << endl;
//  R_out << elem_div(efspnum_matetime,emspnum_old_matetime) << endl;
  R_out << "$Mating.time.effective.Female.new.shell.Spawning.biomass" << endl;
  R_out <<efspbio_new_matetime << endl;
  R_out << "$Mating.time.Female.new.shell.Spawning.numbers" << endl;
  R_out << fspnum_new_matetime << endl;
//  R_out << "$ratio Mating time Female new shell Spawning numbers to male old shell mature numbers " << endl;
//            for(i=styr;i<=endyr;i++){if(emspnum_old_matetime(i)<0.001) emspnum_old_matetime(i)=1.0; }
//  R_out << elem_div(fspnum_new_matetime,emspnum_old_matetime) << endl;
  R_out << "$Predicted.Female.survey.Biomass" << endl;
  R_out << pred_srv1_bioms(1) << endl;
  R_out << "$Predicted.Male.survey.Biomass" << endl;
  R_out << pred_srv1_bioms(2)<< endl;
  R_out << "$Predicted.Female.survey.mature.Biomass" << endl;
  R_out << fspbio_srv1 << endl;
  R_out << "$Predicted.Male.survey.mature.Biomass" << endl;
  R_out << mspbio_srv1<< endl;
  R_out << "$Predicted.total.survey.mature.Biomass" << endl;
  R_out << fspbio_srv1+mspbio_srv1<< endl;
  R_out << "$Predicted.Female.survey.new.mature.numbers" << endl;
  R_out << fspbio_srv1_num(1) << endl;
  R_out << "$Predicted.Female.survey.old.mature.numbers" << endl;
  R_out << fspbio_srv1_num(2) << endl;
  R_out << "$Predicted.Male.survey.new.mature.numbers" << endl;
  R_out << mspbio_srv1_num(1)<< endl;
  R_out << "$Predicted.Male.survey.old.mature.numbers" << endl;
  R_out << mspbio_srv1_num(2)<< endl;

  R_out << "$Observed.Prop.fishery.ret.new.males"<< endl;
  for (i=1; i<=nobs_fish; i++) R_out << yrs_fish(i) << " " << obs_p_fish_ret(1,i)<< endl;
  R_out << "$Predicted.length.prop.fishery.ret.new.males" << endl;
  for (i=1; i<=nobs_fish; i++) 
   {
    ii=yrs_fish(i);  
    R_out <<  ii  <<  " "  <<  pred_p_fish_fit(1,ii)  << endl;
   }
  R_out << "$Observed.Prop.fishery.ret.old.males"<< endl;
  for (i=1; i<=nobs_fish; i++) R_out << yrs_fish(i) << " " << obs_p_fish_ret(2,i)<< endl;
  R_out << "$Predicted.length.prop.fishery.ret.old.males" << endl;
  for (i=1; i<=nobs_fish; i++)
   {
    ii=yrs_fish(i);  
    R_out <<  ii  <<  " "  <<  pred_p_fish_fit(2,ii)  << endl;
   }

  R_out << "$Observed.Prop.fishery.total.new.males"<< endl;
  for (i=1; i<=nobs_fish_discm; i++) R_out << yrs_fish_discm(i) << " " << obs_p_fish_tot(1,i) << endl;
  R_out << "$Predicted.length.prop.fishery.total.new.males" << endl;
  for (i=1; i<=nobs_fish_discm; i++)
   {
    ii=yrs_fish_discm(i);  
    R_out <<  ii  <<  " "  <<  pred_p_fish(1,ii)  << endl;
   }
  R_out << "$Observed.Prop.fishery.total.old.males"<< endl;
  for (i=1; i<=nobs_fish_discm; i++) R_out << yrs_fish_discm(i) << " " << obs_p_fish_tot(2,i) << endl;
  R_out << "$Predicted.length.prop.fishery.total.old.males" << endl;
  for (i=1; i<=nobs_fish_discm; i++)
   {
    ii=yrs_fish_discm(i);  
    R_out <<  ii  <<  " "  <<  pred_p_fish(2,ii)  << endl;
   }
  R_out << "$Observed.Prop.fishery.discard.new.males"<< endl;
  for (i=1; i<=nobs_fish_discm; i++) R_out << yrs_fish_discm(i) << " " << obs_p_fish_discm(1,i) << endl;
  R_out << "$Observed.Prop.fishery.discard.old.males"<< endl;
  for (i=1; i<=nobs_fish_discm; i++) R_out << yrs_fish_discm(i) << " " << obs_p_fish_discm(2,i)<< endl;

  R_out << "$Observed.length.prop.fishery.discard.all.females" << endl;
  for (i=1; i<=nobs_fish_discf; i++) R_out <<  yrs_fish_discf(i)  <<  " "  <<  obs_p_fish_discf(i)  << endl;
  R_out << "$Predicted.length.prop.fishery.discard.all.females" << endl;
  for (i=1; i<=nobs_fish_discf; i++)
   {
    ii=yrs_fish_discf(i);  
    R_out <<  ii  <<  " "  <<  pred_p_fish_discf(ii)  << endl;
   }
  R_out << "$Observed.length.prop.snow.fishery.females" << endl;
  for (i=1; i<=nobs_snowfish_discf; i++)
   {
    R_out <<  yrs_snowfish_discf(i)  <<  " "  <<  obs_p_snow(1,i)  << endl;
   }
  R_out << "$Predicted.length.prop.snow.fishery.females" << endl;
  for (i=1; i<=nobs_snowfish_discf; i++)
   {
    ii=yrs_snowfish_discf(i);  
    R_out <<  ii  <<  " "  <<  pred_p_snow(1,ii)  << endl;
   }
  R_out << "$Observed.length.prop.snow.fishery.males" << endl;
  for (i=1; i<=nobs_snowfish_discm; i++)
   {
    R_out <<  yrs_snowfish_discm(i)  <<  " "  <<  obs_p_snow(2,i)  << endl;
   }
  R_out << "$Predicted.length.prop.snow.fishery.males" << endl;
  for (i=1; i<=nobs_snowfish_discm; i++)
   {
    ii=yrs_snowfish_discm(i);  
    R_out <<  ii  <<  " "  <<  pred_p_snow(2,ii)  << endl;
   }
  R_out << "$Observed.length.prop.redk.fishery.females" << endl;
  for (i=1; i<=nobs_rkfish_discf; i++)
   {
    R_out <<  yrs_rkfish_discf(i)  <<  " "  <<  obs_p_rk(1,i)  << endl;
   }
  R_out << "$Predicted.length.prop.redk.fishery.females" << endl;
  for (i=1; i<=nobs_rkfish_discf; i++)
   {
    ii=yrs_rkfish_discf(i);  
    R_out <<  ii  <<  " "  <<  pred_p_rk(1,ii)  << endl;
   }
  R_out << "$Observed.length.prop.redk.fishery.males" << endl;
  for (i=1; i<=nobs_rkfish_discm; i++)
   {
    R_out <<  yrs_rkfish_discm(i)  <<  " "  <<  obs_p_rk(2,i)  << endl;
   }
  R_out << "$Predicted.length.prop.redk.fishery.males" << endl;
  for (i=1; i<=nobs_rkfish_discm; i++)
   {
    ii=yrs_rkfish_discm(i);  
    R_out <<  ii  <<  " "  <<  pred_p_rk(2,ii)  << endl;
   }

  R_out << "$Predicted.length.prop.trawl.females" << endl;
  for (i=1; i<=nobs_trawl; i++)
   {
    ii=yrs_trawl(i);  
    R_out <<  ii  <<  " "  <<  pred_p_trawl(1,ii)  << endl;
   }
  R_out << "$Observed.length.prop.trawl.females" << endl;
  for (i=1; i<=nobs_trawl; i++) R_out <<  yrs_trawl(i)  <<  " "  <<  obs_p_trawl(1,i)  << endl;
  R_out << "$Predicted.length.prop.trawl.males" << endl;
  for (i=1; i<=nobs_trawl; i++)
   {
    ii=yrs_trawl(i);  
    R_out <<  ii  <<  " "  <<  pred_p_trawl(2,ii)  << endl;
   }
  R_out << "$Observed.length.prop.trawl.males" << endl;
  for (i=1; i<=nobs_trawl; i++) R_out <<  yrs_trawl(i)  <<  " "  <<  obs_p_trawl(2,i)  << endl;

  R_out << "$Observed.Length.Prop.survey.immature.new.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    R_out << ii <<" " <<obs_p_srv1_len(1,1,1,i) << endl;
   }
  R_out << "$Predicted.length.prop.survey.immature.new.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    R_out << ii << " " << pred_p_srv1_len_new(1,1,ii) << endl;
   }
  R_out << "$Observed.Length.Prop.survey.immature.old.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    R_out << ii <<" " <<obs_p_srv1_len(1,2,1,i) << endl;
   }
  R_out << "$Predicted.length.prop.survey.immature.old.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    R_out << ii << " " << pred_p_srv1_len_old(1,1,ii) << endl;
   }
 
  R_out << "$Observed.Length.Prop.survey.immature.new.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++) R_out << yrs_srv1_length(i) <<" " <<obs_p_srv1_len(1,1,2,i) << endl;
  R_out << "$Predicted.length.prop.survey.immature.new.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    R_out << ii << " " << pred_p_srv1_len_new(1,2,ii) << endl;
   }
  R_out << "$Observed.Length.Prop.survey.immature.old.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++) R_out << yrs_srv1_length(i) <<" " <<obs_p_srv1_len(1,2,2,i) << endl;
  R_out << "$Predicted.length.prop.survey.immature.old.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
  {
   ii=yrs_srv1_length(i);  
   R_out << ii << " " << pred_p_srv1_len_old(1,2,ii) << endl;
  }
  R_out << "$Observed.Length.Prop.survey.mature.new.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    R_out << ii <<" " <<obs_p_srv1_len(2,1,1,i) << endl;
   }
  R_out << "$Predicted.length.prop.survey.mature.new.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    R_out << ii << " " << pred_p_srv1_len_new(2,1,ii) << endl;
   }
  R_out << "$Observed.Length.Prop.survey.mature.old.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    R_out << ii <<" " <<obs_p_srv1_len(2,2,1,i) << endl;
   }
  R_out << "$Predicted.length.prop.survey.mature.old.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    R_out << ii << " " << pred_p_srv1_len_old(2,1,ii) << endl;
   }
 
  R_out << "$Observed.Length.Prop.survey.mature.new.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++) R_out << yrs_srv1_length(i) <<" " <<obs_p_srv1_len(2,1,2,i) << endl;
  R_out << "$Predicted.length.prop.survey.mature.new.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    R_out << ii << " " << pred_p_srv1_len_new(2,2,ii) << endl;
   }
  R_out << "$Observed.Length.Prop.survey.mature.old.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++) R_out << yrs_srv1_length(i) <<" " <<obs_p_srv1_len(2,2,2,i) << endl;
  R_out << "$Predicted.length.prop.survey.mature.old.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    R_out << ii << " " << pred_p_srv1_len_old(2,2,ii) << endl;
   }
//for females don't have length data in first four years first year is 1974
     R_out << "$Observed.Length.Prop.survey.all.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    R_out << ii <<" " <<obs_p_srv1_len(1,1,1,i)+obs_p_srv1_len(2,1,1,i)+obs_p_srv1_len(2,2,1,i)<< endl;
              tmpp4+=obs_p_srv1_len(1,1,1,i)+obs_p_srv1_len(2,1,1,i)+obs_p_srv1_len(2,2,1,i);
   }
  R_out << "$Predicted.length.prop.survey.all.females" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    R_out << ii << " " << pred_p_srv1_len_new(1,1,ii)+pred_p_srv1_len_new(2,1,ii)+pred_p_srv1_len_old(2,1,ii) << endl;
    tmpp1+=pred_p_srv1_len_new(1,1,ii)+pred_p_srv1_len_new(2,1,ii)+pred_p_srv1_len_old(2,1,ii);
   }
  R_out << "$Observed.Length.Prop.survey.all.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);
    R_out << ii <<" " <<obs_p_srv1_len(1,1,2,i)+obs_p_srv1_len(1,2,2,i)+obs_p_srv1_len(2,1,2,i)+obs_p_srv1_len(2,2,2,i)<< endl;
         tmpp2+=obs_p_srv1_len(1,1,2,i)+obs_p_srv1_len(1,2,2,i)+obs_p_srv1_len(2,1,2,i)+obs_p_srv1_len(2,2,2,i);
   }
  R_out << "$Predicted.length.prop.survey.all.males" << endl;
  for (i=1; i<=nobs_srv1_length; i++)
   {
    ii=yrs_srv1_length(i);  
    R_out << ii << " " << pred_p_srv1_len_new(1,2,ii)+pred_p_srv1_len_new(2,2,ii)+pred_p_srv1_len_old(2,2,ii) << endl;
  tmpp3+=pred_p_srv1_len_new(1,2,ii)+pred_p_srv1_len_new(2,2,ii)+pred_p_srv1_len_old(2,2,ii);
   }
  R_out << "$Sum.of.predicted.prop.survey.all.females" << endl;
            R_out <<tmpp1<<endl;
  R_out << "$Sum.of.predicted.prop.survey.all.males" << endl;
            R_out <<tmpp3<<endl;
  R_out << "$Sum.of.Observed.prop.survey.all.females" << endl;
            R_out <<tmpp4<<endl;
  R_out << "$Sum.of.Observed.prop.survey.all.males" << endl;
            R_out <<tmpp2<<endl;

  R_out << "$Predicted.mean.postmolt.length.females"<< endl;
  R_out << mean_length(1) << endl;
  R_out << "$Predicted.mean.postmolt.length.males"<< endl;
  R_out << mean_length(2)<<endl; 
  R_out << "$af1" << endl;
  R_out << af1 << endl;
//  R_out << "$af2" << endl;
//  R_out << af2 << endl;
  R_out << "$am1" << endl;
  R_out << am1 << endl;
//  R_out << "$am2" << endl;
//  R_out << am2 << endl;
  R_out << "$bf1" << endl;
  R_out << bf1 << endl;
//  R_out << "$bf2" << endl;
//  R_out << bf2 << endl;
  R_out << "$bm1" << endl;
  R_out << bm1 << endl;
//  R_out << "$bm2" << endl;
//  R_out << bm2 << endl;
  R_out<<"$Predicted.probability.of.maturing.females"<<endl;
  R_out<<maturity_est(1)<<endl;
  R_out<<"$Predicted.probability.of.maturing.males"<<endl;
  R_out<<maturity_est(2)<<endl;
  R_out<<"$molting.probs.female"<<endl;
  R_out<<moltp(1)<<endl;
  R_out<<"$molting.probs.male"<< endl;
  R_out<<moltp(2)<<endl;
  R_out <<"$Molting.probability.mature.males"<< endl;
  R_out <<moltp_mat(2)<<endl;
  R_out << "$observed.pot.fishery.cpue.1979.fishery.to.endyr.fishery" << endl;
  R_out <<cpue(1969,endyr-1)<<endl;
  R_out << "$predicted.pot.fishery.cpue.1978.to.endyr1.survey" << endl;
  R_out <<cpue_pred(styr,endyr-1)<<endl;
  R_out << "$observed.retained.catch.biomass" << endl;
  R_out << catch_ret(1965,endyr-1) << endl;
  R_out << "$predicted.retained.catch.biomass" << endl;
  R_out << pred_catch_ret(styr,endyr-1)<<endl;
  R_out << "$predicted.retained.new.catch.biomass" << endl;
  R_out << (catch_male_ret_new*wtm)(styr,endyr-1)<<endl;
  R_out << "$predicted.retained.old.catch.biomass" << endl;
  R_out << (catch_male_ret_old*wtm)(styr,endyr-1)<<endl;
  R_out << "$observed.retained.discard.male.catch.biomass" << endl;
  R_out << obs_catchtot_biom(1992,endyr-1) << endl;
  R_out << "$predicted.retained.discard.male.catch.biomass" << endl;
  R_out << pred_catch(styr,endyr-1) << endl;
  R_out << "$predicted.retained.discard.new.male.catch.biomass" << endl;
  R_out << (catch_lmale_new*wtm)(styr,endyr-1) << endl;
  R_out << "$predicted.retained.discard.old.male.catch.biomass" << endl;
  R_out << (catch_lmale_old*wtm)(styr,endyr-1) << endl;
  R_out << "$observed.discard.male.mortality.biomass"<<endl;
  R_out << (obs_catchtot_biom(1992,endyr-1)-catch_ret(1992,endyr-1)) <<endl;
  R_out << "$predicted.discard.male.catch.biomass" << endl;
  R_out << pred_catch(styr,endyr-1) -pred_catch_ret(styr,endyr-1)<< endl;
  R_out << "$observed.female.discard.mortality.biomass" << endl;
  R_out << obs_catchdf_biom(1992,endyr-1) << endl;
  R_out << "$predicted.female.discard.mortality.biomass" << endl;
  R_out << pred_catch_disc(1)(styr,endyr-1) << endl;
  R_out << "$observed.male.discard.mortality.biomass" << endl;
  R_out << obs_catchdm_biom(1992,endyr-1) << endl;
  R_out << "$observed.trawl.catch.biomass"<<endl;
  R_out << obs_catcht_biom(yrs_trawl_c)<<endl;
  R_out << "$predicted.trawl.catch.biomass"<<endl;
  R_out <<pred_catch_trawl<<endl;
  R_out << "$observed.snow.female.discard.mortality.biomass" << endl;
  R_out << catch_snowodisc(1) << endl;
  R_out << "$predicted.snow.female.discard.mortality.biomass" << endl;
  R_out << pred_catch_female_snowd << endl;
  R_out << "$observed.snow.male.discard.mortality.biomass" << endl;
  R_out << catch_snowodisc(2) << endl;
  R_out << "$predicted.snow.male.discard.mortality.biomass" << endl;
  R_out << pred_catch_snowd << endl;  
  R_out << "$observed.redk.female.discard.mortality.biomass" << endl;
  R_out << catch_rkodisc(1) << endl;
  R_out << "$predicted.redk.female.discard.mortality.biomass" << endl;
  R_out << pred_catch_female_rkd << endl;
  R_out << "$observed.redk.male.discard.mortality.biomass" << endl;
  R_out << catch_rkodisc(2) << endl;
  R_out << "$predicted.redk.male.discard.mortality.biomass" << endl;
  R_out << pred_catch_rkd << endl;
  R_out << "$predicted.total.male.catch.biomass" << endl;
  R_out<<pred_catch(styr,endyr-1)+pred_catch_rkd(styr,endyr-1)+pred_catch_snowd(styr,endyr-1)+pred_catch_trawl(styr,endyr-1)/2.0<<endl;
  R_out << "$predicted.total.female.catch.biomass" << endl;
  R_out<<pred_catch_disc(1)(styr,endyr-1)+pred_catch_female_rkd(styr,endyr-1)+pred_catch_female_snowd(styr,endyr-1)+pred_catch_trawl(styr,endyr-1)/2.0<<endl;
//  R_out<<"$Estimated total catch div. by male spawing biomass at fishtime"<<endl;
//  R_out<<elem_div(pred_catch(styr,endyr-1)+pred_catch_rkd(styr,endyr-1)+pred_catch_snowd(styr,endyr-1)+pred_catch_trawl(styr,endyr-1)/2.0,mspbio_fishtime(styr,endyr-1))<<endl;
//  R_out << "$estimated retained catch div. by male spawning biomass at fishtime" << endl;
//  R_out <<elem_div(pred_catch_ret,mspbio_fishtime)(styr,endyr-1) << endl;
//  R_out << "$estimated total catch div. by male spawning biomass at fishtime" << endl;
//  R_out <<elem_div(pred_catch,mspbio_fishtime)(styr,endyr-1) << endl;
//  R_out << "$estimated total catch of legal males by legal males at fishtime" << endl;
//   R_out <<elem_div(pred_catch_gt101(styr,endyr-1),bio_males_gt101(styr,endyr-1)) << endl;
//  R_out << "$estimated total catch numbers of males >101 div. by males numbers >101 at fishtime" << endl;
//  R_out <<elem_div(pred_catch_no_gt101(styr,endyr-1),num_males_gt101(styr,endyr-1)) << endl;
//  R_out << "$estimated total catch numbers of males 101 div. by survey estimate males numbers >101 at fishtime" << endl;
//  for(i=styr;i<endyr;i++) obs_tmp(i) = obs_lmales(i-(styr-1));
//  R_out <<elem_div(pred_catch_no_gt101(styr,endyr-1),obs_tmp(styr,endyr-1)*mfexp(-M_matn(2)*(7/12))) << endl;
//  for(i=styr;i<endyr;i++) obs_tmp(i) = obs_lmales_bio(i-(styr-1));

//  R_out << "$estimated total catch biomass of males >101 div. by survey estimate male biomass >101 at fishtime" << endl;
//  R_out <<elem_div(pred_catch_gt101(styr,endyr-1),obs_tmp(styr,endyr-1)*mfexp(-M_matn(2)*(7/12)) ) << endl;

//  R_out << "$estimated total catch biomass div. by survey estimate male mature biomass at fishtime" << endl;
//  R_out <<elem_div(pred_catch(styr,endyr-1),((obs_srv1_spbiom(2))(styr,endyr-1))*mfexp(-M_matn(2)*(7/12))) << endl;

  R_out << "$estimated.annual.total.directed.fishing.mortality" << endl;
  R_out << fmort(styr,endyr-1) << endl;
  R_out << "$estimated.annual.snow.fishing.mortality" << endl;
  R_out << fmortd_snow(styr,endyr-1) << endl;
  R_out << "$estimated.annual.red.king.fishing.mortality" << endl;
  R_out << fmortd_rk(styr,endyr-1) << endl;
  R_out << "$estimated.annual.total.fishing.mortality" << endl;
  for(i=styr;i<=(endyr-1);i++) R_out << F(1,i)(nlenm)+ Fdisct(2,i)(nlenm)+Fdisc_snow(2,i)(nlenm)+Fdisc_rk(2,i)(nlenm) <<" "; R_out<< endl;
  R_out <<"$retained.F" << endl;
  for(i=styr;i<=(endyr-1);i++) R_out <<F_ret(1,i)(nlenm)<<" "; R_out<<endl;
  R_out <<"$ghl" << endl;
  R_out <<catch_ghl/2.2<<endl;
  R_out << "$estimated.annual.fishing.mortality.females.pot" << endl;
  for(i=styr;i<=(endyr-1);i++) R_out << Fdiscf(i)(nlenm) <<" "; R_out<<endl;
  R_out << "$estimated.annual.fishing.mortality.trawl.bycatch" << endl;
  R_out << fmortt(styr,endyr-1) <<endl;
//recruits in the model are 1978 to 2004, the 1978 recruits are those that enter the population
//in spring of 1979, before the 1979 survey - since using the survey as the start of the year
// in the model spring 1979 is stil 1978.  the last recruits are 2003 that come in spring 2004
  R_out << "$estimated.number.of.recruitments.female" << endl;
  for(i=styr; i<=1973; i++) R_out << mfexp(mean_log_rec1_early+rec_devf_early(i))<<" ";
  for(i=1974; i<=(endyr-1); i++) R_out << mfexp(mean_log_rec(1)+rec_dev(1,i))<<" ";
  R_out <<endl<< "$estimated.number.of.recruitments.male" << endl;
  for(i=styr; i<=1973; i++) R_out << mfexp(mean_log_rec1_early+rec_devf_early(i))<<" ";
  for(i=1974; i<=(endyr-1); i++) R_out << mfexp(mean_log_rec(2)+rec_dev(2,i))<<" ";
  for(i=1;i<=median_rec_yrs;i++)R_out<<2*median_rec<<" "; R_out<<endl;
  
  R_out<<"$distribution.of.recruits.to.length.bins"<<endl;
  R_out<<rec_len<<endl;
  R_out<<"$len_bins"<<endl;
  R_out<< "27.5 32.5 37.5 42.5 47.5 52.5 57.5 62.5 67.5 72.5 77.5 82.5 87.5 92.5 97.5 102.5 107.5 112.5 117.5 122.5 127.5 132.5 137.5 142.5 147.5 152.5 157.5 162.5 167.5 172.5 177.5 182.5"<<  endl;
  R_out << "$sel_fsh_m_n"<< endl;
  R_out << sel(1) << endl;
  R_out << "$sel_fsh_m_o"<< endl;
  R_out << sel(2) << endl;
  R_out << "$selectivity.fishery.ret.new.males"<< endl;
  R_out << sel_fit(1) << endl;
  R_out << "$selectivity.fishery.ret.old.males"<< endl;
  R_out << sel_fit(2) << endl;
  R_out <<"$retention.curve.males.new"<< endl;
  R_out <<sel_ret(1,endyr-1)<<endl;
  R_out <<"$retention.curve.males.old"<< endl;
  R_out <<sel_ret(2,endyr-1)<<endl;
  R_out << "$selectivity.discard.females"<< endl;
  R_out <<sel_discf<<endl;
  R_out << "$selectivity.trawl.females"<< endl;
  R_out <<sel_trawl(1,1)<<endl;
  R_out <<sel_trawl(2,1)<<endl;
  R_out <<sel_trawl(3,1)<<endl;  
  R_out << "$selectivity.trawl.males"<< endl;
//  R_out <<sel_trawl(2)<<endl;
  R_out <<sel_trawl(1,2)<<endl;
  R_out <<sel_trawl(2,2)<<endl;
  R_out <<sel_trawl(3,2)<<endl;
  R_out << "$selectivity.snow.females"<< endl;
  R_out <<sel_disc_snow(1,1)<<endl;
  R_out <<sel_disc_snow(2,1)<<endl;
  R_out <<sel_disc_snow(3,1)<<endl;  
  R_out << "$selectivity.snow.males"<< endl;
  R_out <<sel_disc_snow(1,2)<<endl;
  R_out <<sel_disc_snow(2,2)<<endl;
  R_out <<sel_disc_snow(3,2)<<endl;
  R_out << "$selectivity.redk.females"<< endl;
  R_out <<sel_disc_rkc(1,1)<<endl;
  R_out <<sel_disc_rkc(2,1)<<endl;
  R_out <<sel_disc_rkc(3,1)<<endl;  
  R_out << "$selectivity.redk.males"<< endl;
  R_out <<sel_disc_rkc(1,2)<<endl;
  R_out <<sel_disc_rkc(2,2)<<endl;
  R_out <<sel_disc_rkc(3,2)<<endl;  
  R_out << "$srv_sel_m"<<endl;
  R_out << sel_srv_m << endl;
  R_out << "$srv_sel_f"<<endl;
  R_out << sel_srv_f << endl;
  R_out.close();
     
  }
// ===============================================================================

GLOBALS_SECTION
 #include <math.h>
 #include <admodel.h>
  #include <time.h>
  ofstream R_out;
 #undef log_input
 #define log_input(object) CheckFile << "# " #object "\n" << object << endl;
  ofstream CheckFile("Check.Out");
  time_t start,finish;
  long hour,minute,second;
  double elapsed_time;
// ===============================================================================

RUNTIME_SECTION
 //one number for each phase, if more phases then uses the last number
  maximum_function_evaluations 500,1000,3000,3000,5000,5000,5000
  convergence_criteria 1,1,.01,.001,1e-5

TOP_OF_MAIN_SECTION
  arrmblsize = 3000000;
  gradient_structure::set_GRADSTACK_BUFFER_SIZE(4000000); // this may be incorrect in
  // the AUTODIF manual.
  gradient_structure::set_CMPDIF_BUFFER_SIZE(150000000);
  gradient_structure::set_NUM_DEPENDENT_VARIABLES(400);
  time(&start);
  R_out.open("R_results.rep");

FINAL_SECTION
 time(&finish); 
 elapsed_time = difftime(finish,start);
 hour = long(elapsed_time)/3600;
 minute = long(elapsed_time)%3600/60;
 second = (long(elapsed_time)%3600)%60;
 cout << endl << endl << "Starting time: " << ctime(&start);
 cout << "Finishing time: " << ctime(&finish);
 cout << "This run took: " << hour << " hours, " << minute << " minutes, " << second << " seconds." << endl << endl;


